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Executive  Summary 

This  project  was  formulated  in  response  to  the  Munitions  Management  Statement  of 
Need  MMSON-07-04,  “Advanced  technologies  for  detection,  discrimination  and 
remediation  of  munitions  and  explosives  of  concern  (MEC):  UXO  technology”.  Several 
years  ago  it  was  detennined  that  discrimination  performance  using  single-component 
electromagnetic  induction  (EMI)  sensors  in  the  field  had  been  uniformly  poor  relative  to 
expectations,  predominantly  due  to  stringent  requirements  on  positional  accuracy  and 
signal-to-noise  ratio  (SNR).  This  lead  to  the  development  of  a  number  of  multi¬ 
transmitter,  multi-receiver  EMI  sensors  that  could  be  deployed  in  a  cued-interrogation 
mode  so  that  data  could  be  collected  without  moving  the  platfonn.  Initials  tests  with  these 
systems  have  shown  excellent  discrimination  potential.  However,  at  many  sites,  we 
believe  that  the  most  cost-effective  discrimination  strategy  will  be  to  deploy  a  time- 
domain  EMI  one-pass  detection  and  discrimination  system  where  positional  and 
orientation  accuracy,  data-density  and  SNR  are  maximized.  In  this  project,  we  conduct  a 
feasibility  analysis  of  such  a  system.  We  focus  on  (i)  estimating  the  position  and 
orientation  accuracies  that  can  be  achieved  in  a  dynamic  survey  mode;  and  (ii)  on 
simulating  the  discrimination  perfonnance  of  a  number  of  candidate  EMI  sensors  when 
subjected  to  different  position  and  orientation  errors. 

For  sensor  positioning  we  use  the  commercial  off  the  shelf,  Novatel  SPAN  system  which 
comprises  a  Global  Positioning  System  (GPS)  receiver  combined  with  a  Honeywell 
HG1700  tactical  grade  Inertial  Motion  Unit  (IMU)  which  has  a  drift  rate  of  1°  hr'1.  The 
data  are  post-processed  in  a  commercially  available  package  called  Inertial  Explorer. 
Two  positional  tests  were  conducted,  the  first  using  a  Rotating  Amusement  Park  (RAP) 
carousel,  and  the  second  using  an  EM61  towed  array. 

The  RAP  was  instrumented  with  a  precision  optical  encoder  and  the  SPAN  system  was 
installed  on  a  tilt-table  at  a  fixed  position  on  the  RAP.  The  known  distance  of  the  SPAN 
system  from  the  axis  of  revolution  of  the  RAP  ride,  the  tilt  table  angle  and  the  optical 
encoder  information  provided  accurate  ground-truth.  The  GPS  and  IMU  data  from  the 
SPAN  were  augmented  by  an  optical  encoder  attached  to  a  bicycle  wheel  that  was 
attached  to  the  RAP  and  which  was  in  constant  contact  with  the  ground.  For  one  run, 
when  the  system  was  moving  at  a  brisk  2  meters  per  second,  the  position  estimated  from 
the  IMU/GPS  combination  had  a  root-mean-square  error  of  between  1.1  and  1.8  cm. 

For  the  second  position  test,  a  laser  prism  was  mounted  on  the  rear  left  comer  of  an 
EM61  towed-array,  about  1.5  m  away  from  the  Novatel  SPAN  system  which  was 
mounted  on  the  tow-bar  in  line  with  the  center  of  the  towed-array.  The  SPAN-predicted 
position  of  the  prism  was  compared  to  the  prism  position  measured  with  a  Trimble 
SCS930  Universal  Total  Station.  Data  were  collected  at  various  speeds  between  0.2  and 
1 .3  meters  per  second,  on  both  smooth  and  rough  ground,  with  the  array  traveling  along, 
straight  or  curved  paths  and  with  ramps  sometimes  present  that  one  side  of  the  cart  had  to 
travel  over.  Amongst  all  the  tests  the  root-mean-squared  (RMS)  horizontal  error  varied 
from  a  minimum  of  0.6  cm  to  a  maximum  of  1.7  cm.  These  represent  upper  bounds  on 
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the  horizontal  positional  error  of  the  SPAN  system.  If  the  error  was  entirely  due  to  the 
SPAN  heading,  the  RMS  heading  error  would  vary  from  a  minimum  of  0.2°  to  a 
maximum  of  0.6°.  RMS  errors  of  the  3-D  positions  (including  elevation)  varied  from  0.8 
to  3.2  cm.  The  runs  with  larger  errors  typically  exhibited  systematic  biases  in  the 
elevation  angles.  This  indicates  that  there  may  1-2  cm  magnitude  systematic  differences 
in  the  SPAN  elevations  along  adjacent  transect  paths. 

Numerical  investigation  of  simulation  performance  occurred  along  two  fronts.  The  first 
was  semi-analytical  and  based  on  estimating  the  expected  errors  in  the  positions  and 
polarizabilities  of  a  buried  metallic  object.  The  second  was  based  on  monte-carlo 
simulations  of  systems  subjected  to  different  amounts  of  sensor,  position  and  orientation 
noise.  Five  different  sensor  systems  were  simulated,  including  two  Geonics  EM63 
equivalent,  one  with  a  single  receiver  coil,  the  other  with  three  orthogonal  receiver  coils. 
The  remaining  three  systems  were  based  on  the  25  transmitter  /  25  receiver  TEMTADS 
system  developed  at  the  Naval  Research  Laboratory.  The  first  was  deployed  in  the 
standard  cued  mode,  with  the  other  two  deployed  in  dynamic  mode,  one  with  a  single 
large  transmitter  around  the  array,  and  the  other  with  four  1  m  by  1  m  transmitters  that 
fire  sequentially.  Both  the  semi-analytical  and  Monte-Carlo  approaches  demonstrated  the 
excellent  performance  of  the  cued-system.  They  also  revealed  that  the  dynamic  variants 
of  that  system  have  significantly  improved  performance  over  an  EM63  equivalent.  In 
particular,  the  multiple  receiver  coils  makes  the  TEMTADS  more  tolerant  of  position  and 
orientation  errors  than  the  EM63.  The  four-transmitter  TEMTADS  with  5  cm  position 
error  and  2°  orientation  error  outperforms  an  EM63  with  1  cm  and  1°  error. 

We  conclude  that  the  SPAN  system  is  able  to  position  a  towed-array  at  better  than  1-2  cm 
positional  error  and  better  than  0.6  0  orientation  error.  Combined  with  the  fact  that  multi¬ 
sensor  systems  such  as  the  TEMTADS  are  more  tolerant  of  position/orientation  error  we 
hypothesize  that,  at  many  sites,  one-pass  detection  and  discrimination  with  a  suitably 
modified  TEMTADS  would  be  feasible.  Furthermore,  depending  on  anomaly  density, 
this  one  pass  survey  would  be  more  cost-effective  than  the  two  pass  mode  required  when 
the  TEMTADS  is  deployed  in  a  cued-mode. 
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1  PROJECT  BACKGROUND  AND  OBJECTIVES 

This  project  was  formulated  in  response  to  the  Munitions  Management  Statement  of 
Need  MMSON-07-04,  “Advanced  technologies  for  detection,  discrimination  and 
remediation  of  munitions  and  explosives  of  concern  (MEC):  UXO  technology”.  The 
prohibitive  costs  of  excavating  all  geophysical  anomalies  are  well  known  and  are  one  of 
the  greatest  impediments  to  efficient  cleanup  of  unexploded  ordnance  (UXO) 
contaminated  lands  at  department  of  defense  (DoD)  sites.  Therefore,  innovative 
discrimination  techniques  are  required.  Several  years  ago  it  was  detennined  that 
discrimination  performance  using  single-component  electromagnetic  induction  (EMI) 
sensors  in  the  field  had  been  unifonnly  poor  relative  to  expectations,  predominantly  due 
to  stringent  requirements  on  positional  accuracy  and  signal-to-noise  ratio  (SNR)  (e.g. 
Bell,  2005).  This  lead  to  the  development  of  a  number  of  multi-transmitter,  multi-receiver 
EMI  sensors  (e.g.  Berkeley  UXO  Discriminator,  Naval  Research  Laboratory  [NRL] 
cued-interrogation  array  or  TEMTADS)  that  could  be  deployed  in  a  cued-interrogation 
mode  so  that  data  could  be  collected  without  moving  the  platform.  Initials  tests  with  these 
systems  have  shown  excellent  discrimination  potential.  However,  the  need  to  deploy  the 
systems  in  a  cued-mode  can  significantly  increase  the  time  and  costs  of  the  geophysical 
survey.  Recent  results  at  the  Fonner  Lowry  Bombing  and  Gunnery  Range  (FLBGR, 
Billings  et  ah,  2007)  and  the  former  Camp  Sibert  (e.g  Billings  et  ah,  2008)  have 
demonstrated  that,  at  certain  sites,  effective  discrimination  is  possible  by  careful 
interpretation  of  single-component  time-domain  EMI  data  collected  in  a  one-pass 
detection  and  discrimination  mode.  The  discrimination  strategy  was  more  effective  when 
the  time-decay  was  measured  over  a  wider  time  range  (180  microseconds  [ps]  to  25 
milliseconds  [ms]).  At  many  sites,  we  believe  that  the  most  cost-effective  discrimination 
strategy  will  be  to  deploy  a  time-domain  EMI  one-pass  detection  and  discrimination 
system  where  positional  and  orientation  accuracy,  data-density  and  SNR  are  maximized. 
In  this  project,  we  conduct  a  feasibility  analysis  of  a  one-pass  detection  and 
discrimination  system. 

1.1  THEORETICAL  BACKGROUND 

EMI  is  generally  considered  to  be  the  most  promising  technology  for  discriminating 
between  UXO  and  non-UXO  items.  In  the  EMI  method,  a  time  varying  field  illuminates  a 
buried,  conductive  target.  Currents  induced  in  the  target  then  produce  a  secondary  field 
that  is  measured  at  the  surface.  EM  data  inversion  involves  using  the  secondary  field 
generated  by  the  target  for  recovery  of  the  position,  orientation,  and  parameters  related  to 
the  target’s  material  properties  and  shape.  In  the  UXO  community,  the  inverse  problem  is 
usually  simplified  by  assuming  that  the  secondary  field  can  be  accurately  approximated 
as  a  dipole. 

In  general,  time-domain  electromagnetic  (TEM)  sensors  use  a  step-off  field  to  illuminate 
a  buried  target.  The  currents  induced  in  the  buried  target  decay  with  time,  generating  a 
decaying  secondary  field  that  is  measured  at  the  surface. 
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The  time-varying  secondary  magnetic  field  B(/)  at  a  location  r  from  the  dipole  m(7)  is 


where  r  =  r  /|r|  is  the  unit-vector  pointing  from  the  dipole  to  the  observation  point,  /  is 

the  3x3  identity  matrix,  /u0  =  4  nx  10'7  H/m  is  the  pennittivity  of  free  space  and  r  =  |r|  is 
the  distance  between  the  center  of  the  object  and  the  observation  point. 

The  dipole  induced  by  the  interaction  of  the  primary  field  B0  and  the  buried  target  is 
given  by 

m(0=—  M(f)-B0  (2) 

Mo 


where  M(t)  is  the  target’s  polarization  tensor.  The  polarization  tensor  governs  the  decay 
characteristics  of  the  buried  target  and  is  a  function  of  the  shape,  size,  and  material 
properties  of  the  target.  The  polarization  tensor  is  written  as: 


Lx{t)  0  0 

0  L2{t)  0 

0  0  L3(t) 


(3) 


where  we  use  the  convention  that  Li{t\)> L2{t\)> L^) ,  so  that  polarization  tensor 
parameters  are  organized  from  largest  to  smallest.  The  polarization  tensor  components 
are  parameterized  such  that  the  target  response  can  be  written  as  a  function  of  a  model 
vector  containing  components  that  are  a  function  of  target  characteristics. 

The  success  of  dipole  model  based  discrimination  algorithms  depends  on  the  accuracy 
with  which  the  polarization  tensor  can  be  estimated.  With  noisy  and  inaccurately 
positioned  data  it  can  be  difficult  to  constrain  the  depth  and  location  of  the  item,  leading 
to  considerable  uncertainty  in  the  magnitude  and  ratio  of  the  components  of  the 
polarization  tensor  (e.g.  Smith  et  al.,  2004).  For  instance,  when  using  EM61  Hand  Held 
(EM61-HH)  data,  Bell  (2005)  found  that  a  SNR  of  30  decibels  (dB)  and  positional 
accuracy  on  the  centimeter  level  were  required  for  successful  discrimination.  Bell  further 
argues  that  the  data  quality  requirements  of  the  EM61-HH  (17  centimeter  [cm]  coils)  are 
indicative  of  those  required  by  other  sensors  such  as  the  regular  1  meter  (m)  by  0.5  m  coil 
EM61.  This  perception  that  such  stringent  data  quality  objectives  are  a  requirement  for 
advanced  discrimination  has  led  many  groups  to  practically  abandon  efforts  to  develop 
one-pass  detection  and  discrimination  systems  and  instead  concentrate  on  the 
development  of  cued  interrogation  platforms  (e.g.  ESTCP  MM-0601). 

Cued  interrogation  is  one  way  to  optimize  data  quality  because  only  a  limited  area  about 
each  anomaly  needs  to  be  surveyed.  As  a  second  data  collection  effort  is  required,  the 
survey  and  processing  costs  for  UXO  remediation  using  cued-interrogation  will  usually 
be  higher  than  for  one-pass  detection  and  discrimination  (D&D).  For  this  latter  mode  of 
surveying,  two  different  sophisticated  towed  array  systems  have  been  developed  over  the 
past  decade:  the  Multi-sensor  Towed  Array  Detection  System  (MTADS)  developed  by 
the  NRL  and  the  Vehicular  Simultaneous  EM  and  Magnetometer  System  (VSEMS) 


2 


September  2008 


MM-1571:  Next  Generation  Data  Collection  Platform 


Sky  Research,  Inc. 


developed  by  Geo-Centers  (now  SAIC).  Both  systems  use  purposely  built  tow-vehicles 
and  trailers  and  both  deploy  arrays  of  Geonics  EM61  metal  detectors  and  total-field 
magnetometers.  The  main  difference  is  that  the  VSEMS  collects  the  EM  and  magnetic 
data  in  one-pass,  while  the  MTADS  uses  two  passes.  Neither  of  these  systems  can 
achieve  the  data  quality  requirements  suggested  by  Bell  (2005)  and  as  he  points  out,  ‘the 
discrimination  performance  using  the  MTADS  EM61  array  in  the  field  has  been 
uniformly  poor  relative  to  expectations  based  on  results  of  inverting  controlled 
measurements  over  ordnance  and  exploded  shell  fragments  ’. 

To  address  the  issue  of  whether  data  quality  requirements  for  discrimination  are  as 
stringent  as  Bell  (2005)  found  for  the  EM61-HH,  we  first  outline  the  four  main  types  of 
features  extracted  from  the  polarization  tensor  that  are  used  for  discrimination: 

1)  Size  of  the  object:  typically  inferred  from  L\ (?)  at  one  or  more  time-channels,  or 
may  involve  some  combination  of  the  three  polarizations; 

2)  “Two-dimensional”  shape  of  the  object:  typically  inferred  from  the  values  of 
L\(t)  and  L2(t )  at  a  particular  time-channel; 

3)  “Three  dimensional”  shape  of  the  object:  typically  inferred  from  the  values  of  all 

three  polarizations.  For  example,  asymmetric  scrap  has  L2(t\)*  Z3(q)  while 

ferrous  items  with  cylindrical  symmetry  have  Lftfy  L2(ti)  =  Lftf  and  aluminium 
items  with  cylindrical  symmetry  have  L\  (q  )  =  L2(t\)>  Z3  (q ) ; 

4)  Time-decay  characteristics  of  the  object:  typically  inferred  from  the  primary 
polarization  L\(t)  using  ratios  of  time-channels  or  power-law  and/or  exponential 
fits  to  the  decay  shape. 

Recent  results  achieved  with  EM61  and  EM63  sensors  deployed  at  FLBGR  and  the 
Former  Camp  Sibert  (Billings  et  al.,  2007  and  2008)  demonstrated  good  discrimination 
performance  in  the  presence  of  positional  errors  on  the  order  of  2  to  5  cm.  The  analysis 
focused  on  features  1)  and  4)  (the  size  and  time-decay  properties)  derived  from  the 
primary  polarization  tensor.  These  features  were  used  to  guide  a  statistical  classifier  and 
rank  anomalies  according  to  UXO  likelihood.  As  only  the  primary  polarization  was  used, 
the  features  were  more  robust  and  tolerant  of  positional  error  and  sensor/geological  noise, 
than  feature  vectors  extracted  from  the  secondary  and  tertiary  polarizations  (the  “two” 
and  “three”  dimensional  features  above).  The  longer  measurement  time  of  the  EM63 
resulted  in  more  effective  discrimination  due  to  the  significantly  different  late-time 
behavior  of  the  UXO  items  compared  to  the  encountered  clutter.  Therefore,  we  suggest 
that  previous  poor  discrimination  perfonnance  was  partly  due  to  the  use  of  unstable 
polarization  tensor  features  (the  secondary  and  tertiary  polarizations)  and  the  lack  of 
time-decay  information  past  1 .2  ms. 

It’s  apparent  that  the  more  accurate  the  position  and  orientation  infonnation,  the  better 
the  discrimination  performance.  However,  we  assert  that,  at  many  sites,  UXO 
discrimination  will  be  achievable  using  one-pass  detection  and  discrimination.  Therefore, 
we  don’t  believe  that  these  systems  should  be  ignored  in  favor  of  static  cued- 
interrogation  platforms. 
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The  following  issues  need  to  be  considered  when  contemplating  using  an  EMI  sensor  in  a 
one-pass  detection  and  discrimination  mode: 

1)  Positional  and  orientation  accuracy:  In  a  dynamic  mode,  discrimination  relies  on 
accurate  relative  position  between  successive  measurements.  We  will  tackle  this 
issue  using  a  Novatel  SPAN  system  comprising  a  real-time  kinematic  global 
positioning  sytem  (RTK  GPS)  and  a  tactical  grade  inertial  motion  unit  (IMU), 
with  the  dynamical  state  of  the  system  detennined  through  a  Kalman  filter. 

2)  Timing  accuracy:  The  accurate  and  stable  synchronization  of  the  different  sensor 
components  is  a  critical  requirement  for  successful  discrimination  based  on 
dynamically  collected  data.  At  a  forward  traverse  rate  of  1  m/s,  a  timing  error  of  1 
ms  translates  to  an  0.1  cm  error  in  the  position. 

3)  Measurement  time  and  data-density:  The  longer  the  time-range  measured,  the 
further  the  array  will  move  during  a  single  measurement,  the  lower  the  density  of 
the  measurements  and  the  lower  the  SNR  at  the  later  time-channels.  Systems  such 
as  the  TEMTADS  have  the  ability  to  measure  the  time-decay  to  25  ms  after 
transmitter  shut-off,  with  a  single  measurement  with  that  time  range  requiring 
about  100  ms  to  collect.  The  Environmental  Security  Technology  Certification 
Program  (ESTCP)  MM-0504  demonstrations  at  FLBGR  and  Camp  Sibert 
demonstrated  that  the  wider  measurement  range  of  the  EM63  (to  25  ms)  resulted 
in  superior  discrimination  performance  to  the  EM61  (to  1  ms).  It’s  not  clear  as  to 
whether  measurements  are  required  all  the  way  out  to  25  ms  and  some 
compromise  of  measuring  out  to,  say,  8  ms,  may  prove  to  be  more  effective 
(measurements  require  32  ms  and  hence  the  data-density  will  be  about  3  times 
greater  than  at  25  ms  decay). 

4)  Signal  to  noise  ratio:  In  cued-mode,  the  data  are  stacked  to  improve  SNR.  In 
dynamic  mode,  the  choice  of  the  number  of  stacks  involves  a  compromise 
between  improving  SNR  and  blurring  of  spatial  detail.  The  lower  the  number  of 
stacks,  the  less  movement  that  occurs  during  each  measurement.  The  greater  the 
current  in  the  transmitter  the  better  the  SNR,  but  it  also  requires  waiting  longer 
before  taking  the  first  measurement.  Movement  during  the  receive  period  also 
introduces  dynamic  motion  noise  which  lowers  the  SNR.  To  combat  that  noise 
source,  an  extra  reference  sensor  could  be  mounted  at  a  higher  elevation  than  the 
detection  sensors.  The  reference  sensor  will  sense  approximately  the  same 
ambient  noise  as  the  other  receivers,  and  due  to  its  increased  height,  will  have  a 
significantly  reduced  response  from  any  items  of  interest. 

5)  Change  in  background  response  from  the  ground:  During  mobile  data  collection, 
the  geometry  of  the  array  relative  to  the  ground  is  constantly  changing.  Emerging 
work  from  project  Strategic  Environmental  Research  and  Development  Program 
(SERDP)  project  MM- 1573  indicates  that  the  resulting  variation  in  the 
background  response  can  be  significant  and  has  the  potential  to  degrade  the 
accuracy  of  any  polarization  tensors  extracted  from  the  data.  With  knowledge  of 
the  array  geometry  relative  to  the  ground,  variations  in  the  background  response 
can  be  modeled  and  removed  from  the  data. 
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6)  Excitation  of  multiple  axes:  To  be  able  to  discriminate  using  feature  vectors  from 
all  three  polarizations,  the  object  must  be  excited  from  a  wide-range  of  different 
angles.  When  surveying  with  the  MTADS  EM61,  array  data  are  collected  in  two 
orthogonal  excitations  as  a  single  direction  does  not  result  in  significant  excitation 
perpendicular  to  the  direction  of  travel.  To  avoid  the  need  for  two  data  collection 
passes,  a  one-pass  detection  and  discrimination  system  must  utilize  an  intelligent 
transmitter  firing  sequence  so  that  all  principal  axes  of  a  buried  object  are  excited. 

The  main  focus  of  this  research  project  is  on  addressing  point  1)  above:  that  is  on 
optimizing  the  position  and  orientation  accuracy  of  the  EMI  sensor.  We  address  this 
challenge  in  the  next  section  through  the  integration  of  a  GPS  with  an  IMU  using  a 
Kalman  filter  framework.  Specifically,  we  rigorously  test  the  positional  accuracy  of  a 
Novatel  SPAN  system  mounted  on  a  rotating  amusement  park  (RAP)  ride.  Then  in  the 
following  section,  we  conduct  a  more  realistic  positional  accuracy  test  with  the  SPAN 
mounted  on  a  towed-array  that  is  driven  over  both  rough  and  smooth  ground.  Next,  we 
use  survey  design  techniques  and  numerical  simulations  to  characterize  the  expected 
performance  of  different  time-domain  EMI  sensors  deployed  in  a  one-pass  detection  and 
discrimination  mode.  In  the  final  section,  we  provide  recommended  modifications  to 
create  an  alternative  dynamic  deployment  mode  of  the  TEMTADS  cued-interrogation 
sensor. 
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2  OPTIMIZATION  OF  POSITION  AND  ORIENTATION 

ACCURACY 

As  outlined  in  the  last  section,  accurate  position  and  orientation  information  are  required 
for  successful  discrimination  using  a  moving  system.  One  of  the  most  effective  ways  to 
obtain  this  information  is  to  combine  a  GPS  (or  a  laser  positioning  system)  with  an  IMU 
(e.g.  Steinhurst  et  ah,  2005;  Billings  et  ah,  2007).  The  GPS  and  IMU  can  both  provide 
information  on  either  or  both  position  and  orientation.  For  the  GPS,  orientation 
information  can  be  obtained  through  the  use  of  multiple  receivers.  For  the  IMU,  position 
and  orientation  can  be  obtained  through  accelerometers  and  gyroscopes.  To  obtain 
position  information  from  the  IMU  requires  a  double  integration  of  the  IMU  acceleration 
data.  Due  to  the  accumulation  of  errors,  IMU  position  errors  will  quickly  exceed  the 
centimeter-level  position  accuracy  specification.  Frequent  updating  of  the  IMU  position 
is  therefore  needed  to  achieve  the  required  accuracy.  GPS  on  the  other  hand,  can  deliver 
excellent  position  accuracy,  but  has  a  problem  called  “cycle  slips”,  which  are  in  essence 
gross  errors  leading  to  a  discontinuity  in  the  trajectory.  The  combination  of  the  two 
measuring  systems,  therefore,  offers  a  number  of  advantages.  In  the  absence  of  cycle 
slips,  the  excellent  positioning  accuracy  of  differential  GPS  can  be  used  to  provide 
frequent  updates  for  the  inertial  system.  The  inertial  sensors’  orientation  information  and 
the  precise  short-term  position  and  velocity  can  be  used  for  GPS  cycle  slip  detection  and 
correction.  The  objective  of  an  effective  processing  strategy  is  to  minimize  position  and 
orientation  errors  by  combining  IMU  and  GPS  measurements  in  a  manner  that  optimizes 
the  attributes  of  each  measurement  system. 

Steinhurst  et  al.  (2005)  implemented  a  combined  GPS  and  IMU  position/orientation 
system  on  the  GEMTADS  towed-array  platfonn.  They  used  three  Trimble  MS750  GPS 
antennas  that  operate  in  a  RTK  moving  base-mode  to  give  platform  positions  at  20  hertz 
(Hz)  and  platform  orientation  at  10  Hz.  A  Crossbow  VG600CA  Vertical  Gyroscope  was 
used  to  provide  platform  pitch  and  roll  at  a  rate  of  up  to  83  Hz.  To  combine  the  GPS 
angles  and  the  IMU  angles,  the  IMU  data  were  high-pass  filtered  to  retain  the  high 
frequency  variation  of  the  platform  angles,  the  GPS  data  were  low-pass  filtered  to  retain 
the  stable  base  in  the  angles,  and  the  two  were  combined.  To  maintain  RTK  status,  the 
Trimble  moving  base  requires  that  all  three  receivers  continuously  track  a  minimum  of 
four  common  satellites.  In  field  trials,  it  was  found  that  the  RTK  status  would  be  lost  for 
several  minutes  after  trees  or  other  obstacles  obscured  the  sky-view  of  one  or  more  of  the 
receivers.  At  the  sites  surveyed  this  RTK  problem  only  occurred  infrequently  and 
Steinhurst  et  al.  (2005)  estimated  that  the  platform  position  was  accurate  to 
approximately  3  cm  and  the  orientation  to  within  1°. 

The  method  developed  in  Steinhurst  et  al.  (2005)  combines  some,  but  not  all,  of  the 
relative  advantages  of  the  GPS  and  IMU.  For  instance,  the  linear  and  angular 
accelerations  from  the  IMU  were  only  used  for  calculating  the  platform  pitch  and  roll  and 
not  for  determining  the  position  of  the  platform.  A  Kalman  Filter  (e.g.  Kalman,  1960) 
provides  a  computational  framework  for  optimally  combining  the  GPS  and  IMU  data.  In 
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the  next  section  we  briefly  describe  the  Kalman  Filter  algorithm  and  then  describe  an 
experiment  conducted  to  determine  the  accuracy  of  the  Novatel  SPAN  system.  This 
comprises  a  GPS  receiver  combined  with  a  Honeywell  HG1700  tactical  grade  IMU  (see 
http://www.novatel.com/Documents/Papers/HG1700  SPAN62.pdf  for  a  list  of  technical 
specifications),  which  has  a  drift  rate  of  1°  hr'1.  The  data  are  post-processed  in  a 
commercially  available  package  called  Inertial  Explorer. 

2.1  KALMAN  FILTER  THEORY 

A  Kalman  Filter  (KF)  is  a  recursive  data  processing  algorithm  that  estimates  the  state  of  a 
noisy  linear  dynamic  system;  i.e.,  a  vector  x  consisting  of  n  variables  that  describe  some 
interesting  properties  of  the  system.  In  our  case,  the  properties  of  interest  are  the  position 
of  the  sensor,  consisting  of  the  ( x ,  y,  z)  coordinates,  and  the  orientation  (p  of  the  sensor. 
In  practice,  a  KF  has  access  to  several  measurement  sources  that  are  corrupted  with  noise. 
In  our  application,  IMU  measurements  are  made  using  3 -axis  accelerometers  and  rate 
gyros,  and  3-dimensional  position  from  a  GPS.  An  additional  measurement  of  sensor 
velocity  is  made  using  the  bicycle  wheel  encoder.  If  the  measurement  noise  sources 
follow  a  Gaussian  distribution,  which  is  a  good  assumption  in  this  case,  then  the  KF 
estimator  is  statistically  optimal  with  respect  to  any  reasonable  measure  for  optimization 
(Grewal  and  Andrews  1993).  The  KF  processes  all  available  measurements  to  estimate 
the  state,  both  accurate  and  inaccurate  measurements.  It  uses  knowledge  of  the  system 
and  sensor  dynamics,  probabilistic  descriptions  of  the  system  and  measurement  noises, 
and  any  available  data  about  the  initial  values  of  the  state. 

A  KF  can  be  called  a  state  estimator  that  works  on  a  prediction-correction  basis.  This 
means  that  it  computes  a  “belief’  in  a  certain  state  estimate  by  first  making  a  prediction 
based  on  the  dynamics  of  the  system  and  later  correcting  this  prediction  using 
measurements  of  the  system.  In  the  case  of  a  sensor  on  a  RAP  (or  a  towed-array  cart  in 
the  field),  there  is  a  time  history  of  the  motion  and  orientation  of  the  sensor  on  the  RAP 
(or  the  cart).  The  KF  first  makes  a  prediction  of  the  sensor  position  and  orientation  based 
on  the  time  history.  Next,  it  corrects  this  estimate  using  measurements,  such  as  relative 
measurements  of  integrated  acceleration  and  angular  velocity,  and  absolute 
measurements  such  as  GPS  position.  The  main  objective  of  the  KF  is  to  optimize  the 
estimate  of  the  true  state  of  the  system.  In  particular,  it  estimates  the  state  and  gives  a 
measure  of  how  certain  it  is  that  the  state  estimate  is  the  true  state.  What  makes 
estimating  the  state  difficult  is  that  the  state  will  change  over  time  and  both  the  state  and 
the  measurements  are  subject  to  noise. 

Formal  derivations  of  KF  theory  are  lengthy,  complicated,  and  available  in  several  texts 
and  journal  articles  (e.g.,  Kalman  1960;  Maybeck  1979;  Grewal  and  Andrews  1993; 
Welch  and  Bishop  2001).  Here  we  present  a  brief  overview  of  KF  theory  in  order  to 
obtain  some  insight  into  how  to  optimize  its  application. 
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The  KF  estimates  the  conditional  probability  of  being  in  state  xa  given  available 
measurements  z/. . .  z/(.  The  probability  of  being  in  state  x/(  given  observations  z/. . .  z/(  is 
called  the  “belief  ”, 

Bel(xk)  =  P(xk\zi,..  .,Zk).  ^ 

This  belief  definition  can  then  be  split  into  the  prior  belief  Bel'  (xk)  and  the  posterior 
belief  Bet  (xk)  using  Bayes'  rule,  the  theorem  of  total  probability,  and  the  Markov 
assumption.  Then, 


Bel~(xk )  =  P(xk\z1, . . . ,  zfc-i) 

=  J  P(xk\xk^1)Bel+(xk^1)dxk^1 

Bel+(xk )  =  P(xk\zi,  ...,zk) 

_  P(zk\xk)Bel~(xk) 

The  prior  belief  is  the  conditional  probability  of  being  at  state  xk  given  all  the 
measurements  z  up  to  step  k.  The  posterior  belief  is  the  conditional  probability  of  being  at 
state  Xk  given  all  the  measurements  z  up  to  and  including  step  k.  In  order  to  compute  the 
beliefs,  expressions  are  needed  for  the  system  model  P(xk  I  xk  i  )  and  the  measurement 
model  P(xk  I  xa). 


The  KF  computes  the  belief  by  first  computing  the  prior  belief  and  then  computing  the 
posterior  belief.  The  computation  of  the  prior  belief  Bel ~  (xk)  can  be  seen  as  a  prediction 
of  the  state  of  the  system  after  one  time  step.  Without  looking  at  new  measurement 
information,  the  prior  belief  estimates  what  the  state  of  the  system  most  likely  is  after  one 
step.  It  uses  a  model  of  the  system  P(xk  I  xa  1)  and  the  belief  in  what  the  state  was  at  the 
last  time  step  Bet  (xk).  Due  to  noise  in  the  system,  there  may  be  a  prediction  error. 

That  is,  the  prediction  may  be  different  from  the  true  state.  The  computation  of  the 
posterior  belief  Bet  (xk)  can  be  seen  as  a  correction  of  the  state  estimate  resulting  from 
the  prediction.  After  the  KF  has  calculated  the  prior  belief,  new  measurement  infonnation 
gives  direct  information  about  the  true  state  of  the  system.  This  measurement  information 
may  be  used  to  correct  the  predicted  state.  For  this  the  KF  uses  a  model  of  the 
measurements  P(xk  I  xk).  This  model  describes  how  likely  it  is  that  given  a  state  xa  a 
sensor  reading  results  in  measurement  zk.  Given  this  model  and  a  real  measurement,  the 
KF  corrects  the  prior  belief  in  the  state  Bel~  (xk).  Thus,  the  KF  predicts  the  true  state 
of  the  system  after  a  time  step  has  passed  and  it  then  corrects  its  prediction  using 
measurements  that  are  made  of  the  true  state.  The  KF  is  therefore  also  called  a  predictor- 
corrector  state  estimator  (Maybeck  1979). 
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2.2  SPAN  POSITION  AND  ORIENTATION  ACCURACY  TESTING 

We  set  out  to  accurately  measure  the  position  and  orientation  accuracy  achievable  in  a 
dynamic  platform.  Position  errors  of  less  than  1  cm  and  orientation  to  0.1  degree  or  less 
were  targeted.  Initially  a  custom  system  featuring  the  following  components  was 
proposed: 

>  On  board  GPS  with  RTK  base  station. 

>  Custom  inertial  measurement  unit  (IMU). 

>  Wheels  custom  fit  with  optical  encoders  to  provide  distance  measurements. 

>  Trimble  laser-based  robotic  total  station  (RTS)  position  measurement  system. 

>  Custom  hardware  data  acquisition  system  (DAS). 

>  Custom  Kalman  filtering  software  to  optimize  the  data  set. 

To  demonstrate  the  feasibility  of  achieving  the  above  accuracies,  a  system  using  an  off 
the  shelf  GPS  and  IMU  from  Novatel  was  used.  The  SPAN  system,  as  it  is  called,  comes 
with  its  own  custom  Kalman  filtering  software,  Inertial  Explorer  (IE),  enabling  a 
significant  cost  savings  to  be  realized  during  this  initial  phase.  IE  is  optimized  for  1  Hz 
position  estimates,  targeting  vehicle  navigation  applications  but  can  be  used  for 
estimation  of  positions  at  up  to  a  100  Hz  rate. 

The  system  in  Figure  1  was  placed  on  a  rotating  amusement  park  (RAP)  ride  at  a  fixed 
distance  from  the  center  of  the  ride.  The  set  up  is  shown  in  Figure  2.  The  RAP  provides 
the  means  to  repeatedly  place  the  detection  system  over  a  locus  of  points  on  a  circle.  One 
of  the  two  encoders  shown  in  Figure  1  was  used  to  measure  the  instantaneous  angle  of  the 
RAP.  This  angle,  along  with  the  known  distance  between  the  center  of  the  RAP  and  the 
detection  system  provides  a  highly  accurate  estimate  of  the  system  position  at  fixed 
points  in  time. 

Figure  1  shows  a  block  diagram  of  the  system  components  used  to  collect  position  data. 
A  laptop  computer  logs  incoming  data  from  both  the  SPAN  system  and  the  DAS.  Both 
the  SPAN  and  DAS  time  tag  their  respective  data  with  highly  accurate  GPS  times.  This 
yields  the  capability  to  highly  accurately  post-process  the  parallel  data  streams.  The 
SPAN  system  time  stamping  is  accurate  to  the  10s  of  nanoseconds,  and  the  DAS  is 
accurate  to  10  microseconds.  The  SPAN  system  is  used  to  measure  GPS  and  IMU  data. 
The  DAS  is  used  to  measure  optical  encoder  data. 
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Figure  1:  Block  diagram  of  SERDP1571  test  system. 
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FIELD  TEST 


Novate!  SPEC  DAS  Center  of  Axis 

IMU  &  GPS  &  Laptop  of  Rotation 


Carousel  Provides  Changing,  but  Repeatable 
Position  Measurements  to  Determine 
System  Accuracy.  Tilt  Table  and  Variable  Rotation 
Speed  Simulate  Changes  in  Towed  Array 
Acceleration  and  Orientation. 

Figure  2:  Rotating  amusement  park  (RAP)  ride  with  position  measurement  system  components. 

A  second  optical  encoder  was  integrated  into  a  bicycle  wheel,  seen  in  the  lower  left  of 
Figure  2.  The  two  encoders  provide  a  very  accurate  measurement  of  the  location  of  the 
SPAN  system  at  any  given  time.  The  next  stage  of  this  project  was  spent  using  IE  to 
predict  SPAN  location  vs.  time  and  comparing  the  predictions  to  the  known  location 
based  on  the  encoder  data. 


In  addition  to  improving  the  accuracy  of  the  SPAN  location  measurement,  the  bicycle 
wheel  measurement  provides  a  proof  of  concept  for  using  such  a  device  on  the  UXO 
measurement  platform  in  the  future. 
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2.3  DESIGN  OF  ENCODER  ASSEMBLIES 

Figure  3  shows  a  detailed  schematic  of  the  center  encoder  assembly.  The  optical 
encoders  purchased  from  US  Digital  require  a  tight  mechanical  tolerance  for  mounting. 
The  Mylar  encoder  wheel  has  1800  lines  used  to  measure  position  and  is  sensed  with  a 
differential  optical  sensor  that  measures  in  quadrature  using  two  detectors  for  a  resolution 
of  7200  positions  for  one  full  rotation.  The  encoder  wheel  must  be  positioned  to  the 
optical  sensor  within  0.026  inches  in  the  radial  direction  and  0.020  inches  in  the  axial 
direction.  The  methods  used  to  mount  the  encoders  to  both  the  center  shaft  assembly  and 
the  bicycle  wheel  must  maintain  these  tolerances  during  normal  operation. 


Figure  3:  Schematic  of  center  shaft  optical  encoder  assembly. 


In  the  case  of  the  center  shaft  assembly,  an  encoder  mounting  scheme  was  designed  to 
maintain  these  tolerances  using  two  shafts  joined  with  a  flexible  coupling.  As  seen  in 
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Figure  3,  detail  B,  the  shaft  used  to  mount  the  encoder  wheel  utilizes  two  flange  bearings 
and  two  locking  shaft  collars  to  limit  the  radial  and  axial  motion  of  this  shaft  between  the 
mount  plates.  This  confines  the  encoder  wheel  to  a  range  of  motion  that  is  within  the 
tolerances  stated  above  to  prevent  the  wheel  from  colliding  with  the  optical  sensor  and 
damaging  the  fragile  encoder  wheel.  A  flexible  shaft  coupling  is  then  used  to  connect  to 
another  shaft  mounted  to  an  X-Y  translation  stage  which  is  clamped  near  the  center  of  the 
carousel. 


DRIVING  GEAR 


FIXED  X-Y  STAGE  TO 
ADJUST  POSITION  OF 
FIXED  SHAFT 


EG 


ENTIRE  CENTER  SHAFT 
ASSEMBLY  ROTATES 
AROUND  FIXED  SHAFT 


FIXED  CENTER 
SHAFT 


ENCODER  WHEEL  MOUNTS  TO  FIXED 
SHAFT,  ENCODER  OPTICAL  UNIT 
ROTATES  WITH  CENTER  ASSEMBLY 

FLEXIBLE  SHAFT  COUPLING 
ALLOWS  FOR 
SLIGHT  MISALIGNMENT 
BETWEEN  CENTERS 
OF  ROTATION 


ROTATING 


Figure  4:  Photograph  of  center  encoder  assembly  mounted  to  carousel. 

Figure  4  is  a  photograph  showing  the  center  shaft  assembly  mounted  to  the  center  of  the 
carousel.  The  X-Y  stage,  flexible  shaft  coupling,  and  center  shaft  are  clamped  to  a 
cylindrical  part  that  is  stationary.  The  driving  gear  is  rotated  by  a  stationary  motor  and 
drives  a  planetary  gear  that  rotates  the  entire  carousel.  Because  there  is  no  geometric 
center  of  rotation  position  indicated  on  the  carousel,  an  iterative  procedure  was  used  to 
position  the  center  encoder  assembly  by  rotating  the  carousel  and  adjusting  the  X-Y  stage 
to  place  the  shaft  at  the  geometric  center  of  rotation.  The  flexible  coupling  allows  for 
two  shafts  to  have  a  slight  misalignment  and  still  count  the  rotations  accurately  with  the 
optical  encoder. 

Figure  5  is  a  schematic  of  the  bicycle  wheel  encoder  assembly.  This  wheel  is  mounted  to 
the  carousel  at  the  outer  edge  and  rolls  on  the  stationary  platform  tracking  the  motion  of 
the  moving  carousel  with  approximately  19.16  rotations  of  the  bike  wheel  for  one 


13 


September  2008 


MM-1571:  Next  Generation  Data  Collection  Platform 


Sky  Research,  Inc. 


carousel  revolution.  A  mountain  bike  suspension  fork  is  used  to  maintain  contact 
between  the  wheel  and  platform.  The  suspension  fork  is  necessary  because  the  vertical 
clearance  between  the  bottom  of  the  rotating  carousel  and  the  stationary  platfonn  varies 
up  to  75  mm  during  one  revolution  of  the  carousel.  A  16-inch  diameter  wheel  is  used  to 
allow  for  clearance  to  mount  the  encoder  assembly  within  the  limits  of  the  fork. 
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Figure  5:  Schematic  of  bike  wheel  encoder  assembly. 


In  the  case  of  the  bicycle  wheel  encoder  assembly,  the  optical  sensor  of  the  encoder  is 
stationary  and  the  encoder  wheel  is  mounted  to  a  hub  shaft  adapter  that  rotates.  Detail  C 
of  Figure  5  shows  the  hub  shaft  adapter  mounted  to  the  rotating  hub  of  the  wheel.  The 
encoder  optical  sensor  mounts  on  a  plate  this  is  mounted  to  the  disc  brake  mount  bosses 
on  the  suspension  fork.  A  centering  fixture  is  used  during  assembly  to  position  all  the 
components  to  the  correct  position.  The  adjustment  screws  are  then  tightened  and  the 
centering  fixture  is  removed.  This  approach  successfully  allowed  the  encoder  wheel  to 
rotate  within  the  mechanical  limits  described  above  and  not  collide  with  the  optical 
sensor. 
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Figure  6:  Bicycle  wheel  assembly  mounted  to  carousel. 

Figure  6  is  a  photograph  of  the  bike  wheel  encoder  assembly  mounted  to  the  carousel.  A 
rigid  support  structure  is  necessary  to  rigidly  couple  the  wheel  assembly  to  the  carousel  to 
accurately  track  the  rotation.  The  mount  structure  must  also  rigidly  resist  the 
compression  force  of  the  air  shock  as  the  suspension  force  is  compressed.  An  adjustment 
mechanism  is  used  to  position  the  angle  of  the  fork  and  set  the  vertical  position  to 
maintain  rolling  contact  as  the  carousel  rotates. 

2.4  RESULTS 

On  December  13th,  2007,  the  system  was  installed  and  operated.  In  addition  to  the 
system  components  discussed,  a  base  station  GPS  also  was  set  up  and  operated.  The  base 
station  is  necessary  to  attain  the  maximum  accuracy  from  the  SPAN  GPS  measurements. 
The  following  sections  describe  the  results. 

2.4.1  Encoder  Measurements 

Figure  7  is  a  diagram  showing  the  relationship  between  the  center  encoder,  the  wheel 
encoder,  and  the  SPAN  GPS  system.  The  encoder  assemblies  were  installed  on  the 
carousel  and  several  data  collection  runs  were  performed.  The  data  was  quality  checked 
to  insure  smooth  movement  of  the  encoders  and  verify  no  slipping  of  the  center  shaft. 
Figure  8  is  a  plot  of  the  raw  encoder  data.  The  blue  trace  is  the  center  encoder  and  the 
red  trace  is  the  wheel  encoder.  In  this  data  set,  the  carousel  is  gradually  accelerating  and 
coming  up  to  full  speed.  As  expected,  several  wheel  encoder  revolutions  occur  for  each 
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center  encoder  revolution  and  the  peaks  representing  one  revolution  get  closer  together  as 
the  speed  increases. 


PREDICTED  DISTANCE 
TRAVELED  BY  GPS 
BASED  ON 
RGPS  MEASURED 


Figure  7:  Diagram  showing  relationship  between  center  encoder  and  wheel  encoder. 

Based  on  the  information  in  Figure  8,  it  is  possible  to  determine  the  number  of  wheel 
revolutions  for  each  carousel  revolution  and  determine  if  this  number  is  constant,  as 
expected  and  desired.  It  is  also  possible  to  determine  the  distance  traveled  by  the  wheel 
for  each  carousel  revolution.  The  top  graph  in  Figure  9  is  a  plot  of  the  number  of  wheel 
revolutions  for  consecutive  carousel  rotations.  Over  the  course  of  25  rotations  the 
average  number  of  wheel  revolutions  per  carousel  revolution  is  19.169. 
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Figure  8.  Plot  of  center  encoder  and  wheel  encoder  as  carousel  rotates. 


Figure  9.  Number  of  carousel  revolutions  versus  number  of  wheel  revolutions  (top)  and  number 
of  carousel  revolutions  versus  distance  traveled  by  wheel  (bottom). 

The  total  distance  traveled  by  the  wheel  for  one  carousel  revolution  can  be  estimated 
using  the  value  for  Rwheei  arm==142.18  inches  and  multiplying  by  2n  to  get  893.39  inches. 
Since  we  know  the  number  of  revolutions  required  for  the  wheel  to  cover  this  distance, 
we  can  solve  for  the  radius  of  the  bicycle  wheel.  This  is  not  easily  measurable  as  the  air- 
filled  wheel  undergoes  some  compression  due  to  the  force  of  the  suspension  fork. 
Furthermore,  this  compression  changes  depending  on  the  force  applied  by  the  fork. 
Using  an  average  of  19.169  revolutions  to  cover  the  893.39  inches,  the  radius  of  the 
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wheel  is  calculated  to  be  14.835  inches.  Using  this  value,  the  distance  covered  by  the 
wheel  for  each  revolution  is  plotted  in  the  lower  graph  of  Figure  9. 

2.4.2  IMU  Measurements 

On  the  SPAN  system,  the  antenna  was  located  directly  above  the  IMU  center  with  a  zero 
offset  in  the  X  and  Y  direction,  and  0.167  m  of  vertical  separation.  Figure  10  is  a 
photograph  of  the  assembly.  The  IMU  and  GPS  antenna  are  mounted  on  a  tilt  table  to 
vary  the  roll  angle  of  the  assembly.  The  tilt  table  is  precisely  varied  by  rotating  a 
micrometer.  During  one  of  the  experimental  runs,  while  the  carousel  was  stationary,  the 
tilt  table  angle  was  changed  0.763°  by  moving  the  micrometer  two  full  revolutions. 
Figure  1 1  is  a  plot  of  the  roll  angle  recorded  by  the  IMU  at  a  5  Hz  sample  rate.  The  IMU 
roll  angle  indicates  a  change  of  0.768°,  extremely  close  to  the  change  induced  by  the  tilt 
table.  Once  the  carousel  starts  rotating  again,  the  roll  angle  varies  as  the  GPS  rotates 
through  each  carousel  revolution.  The  predicted  vs.  measured  discrepancy  of  0.005° 
highlights  the  accuracy  of  the  SPAN  IMU. 

GPS 

ANTENNA 


IMU 


TILT  TABLE 


Figure  10:  IMU/  GPS  assembly  mounted  to  carousel. 
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Figure  11:  Roll  angle  recorded  by  the  IMU  showing  the  effect  of  varying  the  tilt  table  angle. 


2.4.3  SPAN  Measurements 

Figure  12  shows  a  plot  of  the  estimated  positions  produced  by  the  IE  software  using  the 
SPAN  GPS  and  IMU  data,  plus  the  base  station  GPS  data.  The  location  of  each  data 
point  was  stored  in  Universal  Transverse  Mercator  (UTM)  coordinates,  along  with  the 
time  of  each  data  point,  to  a  data  file.  This  data  was  read  into  Matlab.  The  set  of  data 
points  is  assumed  to  match  the  equation  of  a  circle  centered  at  coordinates  (x0,  y0): 

(x  -  x0)~  +  (y  -  y0)2  =  R2  (6) 

A  least  squares  regression  was  used  to  determine  the  values  of  xG  and  y0,  and  a  new  set  of 
variables  were  solved: 

x'  =  X  -  x0  (7) 

y’=y-Yo  (8) 

Using  these  variables  to  estimate  the  position  on  a  circle  centered  at  x’,  y’  =  0  allowed 
comparisons  between  the  IE  estimated  position  and  the  known  position,  the  latter  based 
on  the  center  encoder  data.  The  least  squares  regression  also  yielded  the  value  of  the 
radius  R:  2.406  meters.  During  the  test  set  up,  the  value  Rgps  was  measured  to  be  2.419 
meters.  This  difference  of  1.33  cm  could  have  been  an  error  in  the  measured  value  of 
Rgps  taken  using  a  tape  measure.  Future  processing  of  more  data  sets,  each  yielding  a 
linear  regression  value  for  the  radius  R,  will  provide  more  confidence  in  the  actual  value. 
However,  for  this  analysis,  the  value  of  2.406  meters  is  used  hereafter. 
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Figure  12:  Plot  of  Inertial  Explorer  estimated  position  output  after  filtering  for  one  of  the 
carousel  measurement  runs. 

Figure  13  (top)  shows  the  estimated  radius  and.  the  estimated  constant  value  versus  time 
during  the  data  collection  period.  For  this  data  set,  the  RAP  was  stopped  for  three 
minutes  during  the  middle  of  the  data  collection,  and  this  period  can  be  seen  to  be 
relatively  quiet  in  the  plot.  On  either  side  of  that  stoppage,  the  carousel  was  moving 
relatively  fast  with  the  SPAN  system  covering  2  meters  each  second. 

Using  the  angle  data  from  the  encoder,  and  the  constant  value  for  the  radius  discussed 
above,  a  known  position  vs.  time  was  generated.  The  absolute  value  of  the  distance  at 
each  point  between  the  known  position  and  the  IE  estimated  position  is  shown  in  the 
bottom  plot  of  Figure  13.  This  is  not  simply  the  difference  between  the  constant  radius 
value  and  estimated  radius  value  shown  in  the  top  plot;  it  also  takes  into  account  the 
angular  position  of  the  SPAN  system  at  each  point  in  time. 

For  the  first  part  of  the  plot  the  mean  error  is  1.54  cm,  the  standard  deviation  is  0.80  cm 
and  the  root-mean-square  error  (RMSE)  (which  combines  the  bias  and  variance)  is  1.74 
cm.  For  the  last  part  of  the  plot  the  mean  is  1.14  cm,  the  standard  deviation  is  0.58  cm 
and  the  RMSE  is  1.28  cm.  The  error  is  less  after  the  stoppage  and  we  suspect  that  this 
may  be  due  to  better  satellite  coverage. 

The  RMSE  in  position  is  between  1.1  and  1.8  cm:  this  is  approaching  our  specified 
performance  objective  of  1  cm.  We  are  encouraged  by  this  result,  as  it  was  obtained  using 
default  values  for  the  noise  densities  of  the  IMU  and  GPS.  Unlike  GPS,  inertial 
computations  are  very  sensitive  to  the  input  of  a  priori  stochastic  information  that  is 
available  for  each  inertial  system.  According  to  Novatel,  as  a  rule  of  thumb,  the  more 
expensive  the  IMU,  the  less  likely  it  is  to  achieve  poor  results  from  a  non-optimal  choice 
of  input  statistical  quantities.  A  poor  choice  of  noise  densities  for  a  navigation  grade 
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inertial  system  may  simply  provide  the  user  with  a  sub-optimal  set  of  outputs,  say  for 
attitude.  An  inappropriate  choice  of  noise  densities  for  an  inexpensive  IMU  (e.g.,  MEMS 
array  with  a  drift  of  300°  per  hour)  may  result  in  meaningless  attitude  information.  The 
Honeywell  tactical  grade  IMU  used  in  the  SPAN  system  is  state-of-the-art  and  has 
excellent  performance  specifications. 

Rather  than  focusing  on  optimizing  the  accuracy  of  the  positions  in  the  RAP  test,  we 
instead  designed  another  more  realistic  test  of  the  positional/orientation  accuracy  of  the 
SPAN. 


Radius 


Figure  13:  Radius  estimate  vs.  time  (top)  and  location  estimate  error  (bottom)  for  one  data  set 
taken  on  the  RAP  platform. 
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3  POSITIONAL  AND  ORIENTATION  ACCURACY  TEST  ON  A 

TOWED  ARRAY 


While  the  accuracies  achieved  with  the  RAP  test  were  encouraging,  the  test  was  not 
representative  of  the  dynamic  environment  of  a  towed-array  or  cart  collecting  data  over  a 
real-world  site.  The  dynamics  in  the  RAP  were  not  representative  of  cart-motions  in  the 
field,  but  perhaps  more  importantly,  the  RAP  tests  concentrated  on  detennining  the 
accuracy  of  predicting  the  position  of  the  GPS  antenna.  In  practice,  there  is  almost  always 
an  offset  (sideways,  forwards,  up/down  or  all  three)  between  the  antenna  and  the  sensor. 
That  is,  for  rectangular  EMI  sensor  we  are  more  interested  in  the  positions  of  the  four 
comers  of  the  transmitter  and  receiver  coils  than  in  the  position  of  the  GPS  antenna.  We 
therefore  decided  to  conduct  a  more  realistic  positional  test  by  mounting  the  SPAN  on 
one  of  Sky  Research’s  EM61  towed  arrays.  “Ground-truth”  information  was  collected  by 
a  laser-based  positioning  system  tracking  a  prism  mounted  on  the  rear-left  comer  of  the 
towed  array  (one  of  the  corners  of  the  transmitter/reciever).  The  laser-prism  has  a  lever- 
arm  of  close  to  1.5  m  relative  to  the  SPAN  system.  Roughly  speaking,  an  error  of  1°  in 
orientation  will  translate  to  approximately  2.5  cm  error  in  the  predicted  location  of  the 
prism.  Therefore,  if,  for  instance,  2.5  cm  accuracy  were  demonstrated,  it  would  place 
upper  bounds  of  2.5  cm  position  and  1°  orientation  error  of  the  SPAN  system.  While  this 
test  does  not  have  the  accuracy  and  rigor  of  the  RAP  test,  it  does  have  the  potential  to 
estimate  the  positional  and  orientation  accuracy  under  realistic  data  collection  conditions. 

For  the  laser-based  positioning  system  we  used  Trimble’s  SPS930  Universal  Total 
Station,  capable  of  collecting  data  at  20  Hz  with  an  accuracy  of  ±5  mm  +2  mm  parts  per 
million  (ppm).  The  SPS930  was  released  in  August  2007  and  fully  characterized  as  part 
of  SERDP  MM- 1441  (Foley  et  al.  2008).  It  uses  Trimble’s  MultiTrack  technology, 
allowing  the  instrument  to  track  a  passive  prism  (mirror-only)  or  active  prism  (mirrors 
with  an  infrared  (IR)  light  source),  and  Trimble’s  MagDrive  fourth  generation  servo 
technology.  The  MagDrive  technology  uses  magnetic  levitation  to  eliminate  friction  in 
the  system  producing  smooth  horizontal  and  vertical  angular  measurements.  The  SPS930 
uses  a  2.4  gigahertz  (GHz),  115,200  baud  radio  to  receive  data  and  transmit  commands. 
For  this  survey,  we  used  an  active  MT900  prism. 

3.1  METHODOLOGY 
3.1.1  Test-plan 

Figures  14  and  15  show  the  Sky  Research  EM61  towed  array  with  the  Novatel  SPAN 
system  mounted  on  the  tow-bar  support  and  the  MT900  active  prism  on  the  rear-left 
comer  of  the  array.  Tests  were  conducted  using  the  following  general  test-plan 
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1)  Accurately  measure  the  position  of  the  MT900  prism  relative  to  the  SPAN.  The  prism 
position  needs  to  be  determined  at  sub-centimeter  accuracy  relative  to  the  SPAN 
location  for  this  test  to  be  meaningful. 

2)  SPAN  to  SPS930  time-alignment:  Lift  the  right-hand  side  of  the  array  up  by  1-2  feet 
while  collecting  data  with  all  instruments.  Then  lift  the  rear  of  the  array  up  also  while 
collecting  data  with  all  instruments.  These  tests  will  allow  the  time-offset  of  the 
SPS930  relative  to  the  SPAN  to  be  detennined. 

3)  “Parking  lot”  tests:  On  a  smooth  surface  perfonn  the  following  tests: 

a.  Drive  a  50  m  long  straight  line  at  approximately  0.5  miles  per  hour  (mph)  then 
repeat  at  1,  1.5,  2  and  3  mph  . 

b.  Drive  in  an  S-pattern  along  the  same  50  m  section  at  0.5,  1,  1.5,  2  and  3  miles 
per  hour. 

c.  Put  several  small  “speed  bumps”  along  the  path  such  that  only  the  left  wheel 
will  go  over  them  and  drive  in  a  straight  line  at  0.5  and  1  mph  (try  and  do  1.5, 
2  and  3  mph  if  the  array  and  gator  can  handle  the  bumps. 

d.  Repeat  c  with  the  speed  bumps  on  the  right-hand  wheel. 

4)  “Rough  field”  tests:  Repeat  the  tests  in  3)  on  a  rough  field. 

In  early  June  2008  we  made  an  initial  attempt  to  conduct  the  tests  described  above  at  the 
FLBGR.  Mobilization  and  equipment  rental  costs  were  minimized  by  conducting  the  test 
in  conjunction  with  a  test  of  the  Sky  Research  HeliMag  system  as  part  of  ESTCP  MM- 
0741.  The  SPAN  data  were  logged  on  an  internal  flash-card  on  the  SPAN  receiver  and 
there  was  no  fonnal  integration  of  the  SPAN  and  DAS.  A  second  GPS  receiver  on  the 
towed-array  was  used  to  provide  the  timing  synchronization  pulse  to  the  DAS.  The  test 
was  essentially  a  total  failure:  the  SPAN  kept  rebooting  due  to  a  faulty  power-connector 
and  the  problem  went  undetected  by  the  operator  as  no  information  on  the  state  of  the 
SPAN  was  displayed  in  the  software  controlling  the  DAS.  The  project  team  had  to 
demobilize  before  the  power-problem  could  be  remedied.  We  decided  to  tightly  integrate 
the  SPAN  into  the  DAS  before  attempting  another  test. 
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Figure  14:  Sky  Research  EM61  towed  array  with  GPS  receiver  (point  81),  HG1700  IMU  (points 
4  and  9)  and  MT900  prism  (point  132).  The  numbers  delineate  reference  points  used  to  define  the 
precise  geometry  of  the  array  through  analysis  of  multiple  photographs  taken  from  different 
positions. 


Figure  15:  The  new  Trimble  SPS930  (left)  uses  radio  modems  to  send  commands  and  positions 
data  to  and  from  a  portable  controller  while  tracking  the  prism  on  the  bottom  comer  of  the  towed 
array. 
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3.1.2  System  integration 

To  prevent  a  repeat  of  the  system  failure  experienced  during  the  first  test  at  FLBGR,  the 

SPAN  was  tightly  integrated  into  the  DAS.  This  system  integration  task  required  the 

following: 

•  Characterize  the  reason  for  frequent  resets  of  the  SPAN:  The  SPAN  would  over  load 
the  serial  port  when  streaming  data  to  the  SPAN’S  COM1  while  using  the  Novatel 
control  and  display  unit  (CDU)  executable.  The  Novatel  CDU  executable  was 
initially  used  in  order  to  configure  and  send  commands  to  the  SPAN.  Hyperterminal 
was  finally  used  instead  to  mitigate  the  overhead  and  resulting  overload  of  the  SPAN 
COM  port  of  using  CDU.  Additionally,  a  SPAN  power  supply  interface  was  loose 
and  caused  intennittent  resetting.  This  interface  was  removed  and  replaced  with  a 
direct  connection  to  the  12V  battery. 

•  Characterize  why  the  current  SPAN  configuration  changed  dynamically  and 
seemingly  randomly:  when  using  the  Novatel  CDU  and  exiting,  the  SPAN  would  be 
left  configured  differently  than  just  prior  to  using  the  Novatel  CDU.  After  switching 
to  Hyperterminal  for  command  and  control  this  problem  was  mitigated. 

•  Characterize  why  the  SPAN  lost  the  stored  configuration  in  non-volatile  memory: 
The  SPAN  would  frequently  loose  its  stored  configuration.  The  root  cause  was  a  bad 
cable  between  the  SPAN  and  the  IMU.  When  this  cable  was  rebuilt  the  SPAN 
maintained  the  configuration  in  non-volatile  memory. 

•  A  sequence  of  SPAN  commands  had  to  be  derived  for  our  deployment  needs 
resulting  in  the  following  script: 

FRESET 

COM  COM1  115200  N  8  1  OFF 
INTERFACEMODE  AUX  IMU  IMU  OFF 
COM  COM2  38400  N  8  1  OFF 
INTERFACEMODE  COM2  CMR  NONE 
SETIMUTYPE  IMU_HG1 7 00_AG1 1 
SETIMUORIENTATION  5 
VEHICLEBODYROTATION  000 
SAVECONFIG 
GROUP  CLEAR 
GROUP  ADD  POWERUP 

GROUPLOG  ADD  POWERUP  COM1  INSPVAB  ONTIME  1 
GROUPLOG  ADD  POWERUP  COM1  BESTPOSB  ONTIME  1 
GROUPLOG  ADD  POWERUP  FILE  RAWEPHEMB  ONCHANGED  0 
GROUPLOG  ADD  POWERUP  FILE  RANGECMPB  ONTIME  .1 
GROUPLOG  ADD  POWERUP  FILE  RAWIMUSB  ONNEW 
GROUPLOG  ADD  POWERUP  FILE  INSPVASB  ONTIME  1 
GROUPLOG  ADD  POWERUP  FILE  BESTPOSB  ONTIME  1 
GROUPUSE  START  POWERUP 

•  A  pulse-per-second  (PPS)  interface  had  to  be  created  to  condition  the  PPS  signal  from 
the  SPAN’S  input/output  (I/O)  port  so  that  it  was  compatible  with  the  DAS. 

•  The  second  PPS  port  on  one  of  our  DAS  units  was  modified  to  be  compatible  with  the 
external  PPS  conditioning  circuit. 
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•  Software  processing  modules  were  written  that  would  remove  the  timestamps  and 
stitch  the  binary  blocks  together  and  maintain  byte  alignment. 

•  Various  SPAN/IMU  states  were  exposed  to  the  operator  by  extracting  values  from  the 
SPAN  binary  stream  (Figure  16). 


Debug  Middle  End  Terminal  Field  Notes 


SPASol°sJiy1t;  pTsS  GPSF  ■  EMI  EM2  EM3  DA|Jf1'1  NoEiiois  AC  pSize  BuHer 


FM1  FM7  FM3  IMI I  r - 1 

Figure  16:  Partial  screen  shot  of  the  DAS  user-interface  showing  the  SPAN  system  status. 

•  In  order  to  collect  quality  data  with  the  SPAN,  requirements  and  constraints  were 
discovered  that  must  be  met: 

o  A  sufficient  number  of  satellite  views  for  a  good  GPS  solution 
o  IMU  data  to  start  streaming  after  a  good  GPS  solution  has  been  met 
o  Visual  validation  that  the  SPAN  COM1  and  AUX  back  panel  lights  are  blinking 
and  the  SPAN  LCD  display  is  stating  “Logging” 
o  An  IMU  alignment  to  succeed 

o  Once  the  alignment  has  succeeded,  start  moving  with  enough  heading  and  height 
changes  for  the  IMU  to  converge,  change  state,  and  report 
INSSOLUTIONGOOD 

o  At  this  time  the  vehicle  needs  to  be  stationary  for  at  least  2  minutes  before  the  test 
is  initiated  and  prior  to  the  start  of  data  collection  with  binary  data  being  written 
to  the  SPAN  Flash  card  (this  is  a  requirement  from  the  post  processing  tool 
“Inertia  Explorer”) 

o  After  the  test  run  has  completed  and  “Close  Survey”  has  been  pressed,  the  vehicle 
will  again  have  to  be  stationary  for  2  minutes. 

•  In  order  to  help  the  operator  validate  and  conform,  the  DAS  Front  End  was  modified 
with  a  custom  EM6 1/SPAN  configuration  selectable  from  the  configuration  pull 
down  menu.  This  included  display  of  the  indicators  above  as  well  as  strip-charts  and 
text-fields  of  the  GPS  and  IMU  data. 


3.1.3  Precise  measurement  of  the  array  geometry 

To  meaningfully  compare  the  SPAN  GPS  positions  with  the  SCS930  reported  positions 
requires  precise  measurement  of  the  geometry  of  the  MT900  prism  relative  to  the  GPS 
antenna  phase  center.  Lever  arm  measurements  were  performed  on  the  Sky  towed  array 
system  by  utilizing  close  range  photogrammetry  techniques.  Photomodeler  Scanner, 
version  6.2.2.596  was  used  to  process  the  imaging  data  and  generate  3  dimensional  (3D) 
coordinates. 

Images  of  the  towed  array  were  taken  with  a  Sony  A-200,  10.2  mega  pixel,  digital  single 
lens  reflex  (DSLR)  camera  with  an  18mm  focal  length  lens.  The  A-200  has  a  field  of 
view  crop  factor  of  1.5,  thus  the  field  of  view  for  the  digital  images  is  comparable  to  a 
27mm  focal  length  lens  if  a  full  frame  camera  format  such  as  35mm  film  was  used.  The 
Sony  A-200  was  calibrated  to  determine  finite  values  for  the  focal  length,  lens  distortion, 
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format  aspect  ratio,  and  principal  point.  The  camera  is  calibrated  by  taking  multiple 
images  of  a  calibration  grid  from  different  perspectives  and  processed  using  a  calibration 
routine  in  the  software.  The  computed  finite  camera  parameters  are  stored  in  a  calibration 
fde  to  be  used  during  photo  processing.  The  grid  used  to  calibrate  the  A-200  camera  is 
shown  in  Figure  17. 


C 


•  *•*  • 


Figure  17:  Grid  used  for  calibrating  the  A-200  camera  in  the  Photomodeler  software. 

Four  images  of  the  towed  array  were  taken  and  imported  into  the  photomodeler  software 
(Figure  18).  The  raw  images  were  then  idealized  using  the  camera  calibration  parameters. 
Idealizing  the  image  removes  the  effects  of  lens  distortion,  non-centered  principal  point 
and  non-square  pixels. 

Once  the  photos  have  been  imported  and  idealized  common  points  among  the  photos  are 
selected  and  used  to  reference  and  orient  the  photos  to  each  other.  Figures  14  and  19 
show  examples  of  selected  reference  points  on  the  GPS,  IMU  and  MT900  prism. 

After  the  photos  have  been  referenced  to  each  other  they  are  processed  through  an  image 
bundle  adjustment.  Once  the  bundle  adjustment  has  been  successfully  computed  the 
camera  stations  are  oriented  to  the  image  planes.  Figure  20  displays  a  picture  of  the 
solved  camera  stations  and  image  planes. 

From  the  computed  camera  stations  and  image  planes  a  3D  model  is  generated.  Rotation, 
translation  and  scale  factor  parameters  are  applied  to  the  3D  model  to  allow  for  the 
extraction  of  physical  distances.  The  scale  factor  is  applied  by  defining  a  distance 
between  two  reference  points  in  the  photos.  The  left  panel  of  Figure  20  is  an  image  for 
the  3D  model  points  and  solved  camera  stations. 
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Figure  18:  These  four  photos  of  the  towed  array  were  imported  into  the  photomodeller  software. 


Figure  19:  Example  reference  points  on  MT900  prism,  GPS  and  IMU. 

The  results  of  the  bundle  adjustment  yielded  sub  millimeter  results,  quality  statistics  of 
the  bundle  adjustment  are  listed  below: 


Quality 

Point  Marking  Residuals 

Photographs 

Overall  RMS:  0.402  pixels 

Total  Number:  4 

Maximum:  1.533  pixels 

Bad  Photos:  0 

Point  124  on  Photo  1 

Weak  Photos:  0 

Minimum:  0.061  pixels 

OK  Photos:  4 

Point  72  on  Photo  4 

Number  Oriented:  4 

Maximum  RMS:  1.258  pixels 

Number  with  inverse  camera  flags  set:  0 

Point  124 

Minimum  RMS:  0.058  pixels 

Point  72 
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Cameras 

Camera  1:  DSLR-A200  [18.00] 

Calibration:  yes 

Number  of  photos  using  camera:  0 

Average  Photo  Point  Coverage:  83% 

Camera2:  DSLR-A200  [18.00]  -  Idealized 

Calibration:  yes 

Number  of  photos  using  camera:  4 

Average  Photo  Point  Coverage:  83% 

Point  Tightness 

Maximum:  2.8  mm 

Point  124 

Minimum:  0.082  mm 

Point  72 

Photo  Coverage 

Number  of  referenced  points  outside  of  the  Camera's  calibrated 
coverage:  0 

Point  Precisions 

Overall  RMS  Vector  Length: 
0.817  mm 

Maximum  Vector  Length:  1.08 
mm 

Point  4 1 

Minimum  Vector  Length:  0.653 
mm 

Point  62 

Maximum  X:  0.761  mm 

Maximum  Y :  0.739  mm 

Maximum  Z:  0.564  mm 

Minimum  X:  0.362  mm 

Minimum  Y:  0.276  mm 

Minimum  Z:  0.302  mm 

Figure  20:  Image  planes  and  camera  station  solution  for  the  towed-array  (left)  along  with  a  3-D 
model  of  the  points  and  camera  stations  (right). 


The  quality  statistics  indicate  sub  millimeter  precisions.  However,  after  independently 
verifying  known  distances  in  the  3D  model  an  accuracy  of  approximately  2mm  is 
expected.  This  discrepancy  can  be  attributed  to  the  distance  that  was  applied  to  the  photo 
reference  points  used  to  compute  the  project  scale  factor.  The  3D  model  was  augmented 
in  AutoCAD  and  center  points  for  the  GPS  antenna,  ATS  prism  and  IMU  center  were 
computed  from  the  reference  points  generated  in  the  photo  modeler  software.  The 
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geometry  of  the  SPAN  IMU  and  MT900  prism  relative  to  the  phase-center  of  the  GPS 
antenna  are  given  in  Table  1. 


Position 

Left-right 
offset  (cm) 

Front-back 
offset  (cm) 

Up-down 
offset  (cm) 

SPAN  IMU 

0 

0 

-16.7 

SCS930  prism 

-148.3 

-40.8 

-51.6 

Table  1:  Position  of  the  IMU  and  SCS930  prism  relative  to  the  GPS  receiver  in  the  local 
coordinate  system  of  the  towed-array. 

3.2  RESULTS 

The  SCS930  laser-gun  was  set-up  over  Control  Point  2  (CP-2,  Table  2).  The  MT900 
prism  was  located  over  points  CP-1  and  CP-3  and  back-site  measurements  were  taken  so 
that  distance  and  angle  measurements  from  the  SCS930  could  be  converted  to  UTM, 
NAD-83  coordinates.  The  SCS930  orientation  was  taken  as  the  average  of  the  computed 
orientations  for  the  CP-1  and  CP-3  back-sites.  Raw  distances  measured  between  the 
SCS930  and  the  MT900  prism  were  corrected  for  distortions  in  the  UTM  grid  using  a 
scale  factor  of  0.999330  (that  is  the  distance  measurements  were  multiplied  by  this 
number).  The  SCS930  was  operated  at  a  rate  of  10  Hz. 

The  SPAN  GPS  corrections  came  from  a  permanent  GPS  base-station  on  control  point 
Sky-Base2,  about  9.0  km  from  the  FLBGR  site  compound  where  the  experiments  were 
conducted.  The  GPS  data  were  logged  at  20  Hz  and  the  IMU  at  100  Hz. 


Control  Point 

Easting  (m) 

Northing  (m) 

Elevation  (m) 

CP-1 

527385.316 

4385648.348 

1800.825 

CP-2 

527449.999 

4385652.390 

1802.377 

CP-3 

527480.617 

4385621.359 

1805.264 

Sky-Base2 

536552.830 

4385850.456 

1814.185 

Table  2:  Local  survey  control  around  the  FLBGR  site-compound  along  with  the  main  GPS  base 
station  position  for  the  FLBGR  site.  Positions  are  given  in  UTM  Zone  13N,  NAD-83  (Conus) 
datum,  with  the  Geiod03  (Conus)  geoid. 


3.2.1  Position  and  time-calibrations 

The  SPAN  positions  are  precisely  synchronized  to  GPS  time.  A  1  PPS  output  from  the 
SPAN  disciplines  the  clock  in  the  DAS  which  allows  the  arrival  of  the  SCS930  to  be 
time-stamped  with  an  accuracy  of  10  ps  (actually  plus  or  minus  1  integer  second 
depending  on  the  relative  time  of  arrival  of  the  pulse  UTC  message).  The  SCS930  reports 
its  position  with  a  very  stable  time-lag  that  appears  to  vary  by  less  than  4  milliseconds. 
To  accurately  compare  the  SPAN  and  SCS  reported  positions  we  need  to  precisely 
determine  the  latency  of  the  SPS  positions  relative  to  GPS  time.  We  lift  the  rear-left 
comer  of  the  array  in  a  calibration  measurement,  and  then  vary  the  time-lag  of  the  SCS 
until  it  best  matches  the  positions  reported  by  the  SPAN  (Figure  21).  This  procedure 
revealed  best  fit  lags  for  the  different  coordinates  of  -1.07  seconds  for  easting,  -1.039  for 
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northing  and  -1.047  for  elevation.  There  was  very  little  movement  in  the  easting  so  we 
use  the  average  of  the  northing  and  elevation  channels  to  specify  a  best-fit  lag  of  -1.043 
seconds.  This  lag  is  used  for  all  the  comparisons  that  follow. 


Figure  21:  Comparison  of  SPAN  and  SCS930  positions  during  a  calibration  measurement  where 
the  left-hand  side  of  the  array  is  lifted.  The  SCS930  times  are  corrected  by  -1.043  seconds  for  the 
comparisons.  This  time-lag  was  obtained  through  cross-correlation  of  the  time-series  of  the 
SPAN  and  SCS930  positions  (bottom  right  plot). 

3.2.2  Smooth  surface  tests 

We  start  by  analyzing  the  data  collected  while  the  towed-array  traversed  a  smooth  gravel 
road.  There  were  four  different  types  of  tests  (Figure  22):  (i)  driving  along  a  straight  line; 
(ii)  driving  in  an  S-pattern;  (iii)  driving  straight  with  the  left-side  of  the  array  going  over 
a  ramp;  and  (iv)  same  as  (iii)  but  with  the  right-side  going  over  the  ramp.  Data  were 
collected  at  speeds  of  0.5,  1.0,  1.5,  2.0  and  3.0  mph  for  each  of  these  four  tests. 

The  SPAN  data  were  post-processed  in  Inertial  Explorer  and  the  position  of  the  MT900 
prism  was  predicted  using  the  offsets  measured  by  the  photogrammetry  survey  (Table  1). 
The  SCS930  position  time-base  was  corrected  by  the  measured  time-lag  of  -1.043 
seconds,  and  the  SCS930  positions  were  interpolated  to  the  same  times  as  the  SPAN.  The 
error  in  the  position  prediction  of  the  SPAN  was  then  calculated  by  subtracting  the 
SCS930  measured  positions.  Immediately,  we  noticed  that  there  was  a  bias  in  the 
coordinate  systems  of  the  SPAN  and  SCS930  reported  positions  (Figure  23a).  The 
easting  coordinate  bias  varied  from  about  -1  cm  at  the  most  westerly  point,  to  +4  cm  at 
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the  most  easterly  point.  The  northing  and  elevation  biases  were  approximately  constant 
with  values  of  +7  cm  and  -1772  cm  respectively  (the  SPAN  and  SCS  height-datums  were 
different).  The  bias  could  result  from  one  or  more  of  the  following  issues:  (i)  Offset 
between  the  SPAN  and  SCS930  control  points  (translation  error);  (ii)  incorrect 
establishment  of  the  orientation  of  the  SCS930  through  the  back-site  procedure  (rotation 
error);  and/or  (iii)  incorrect  scale  factor  (dilation  error).  To  correct  for  this  coordinate 
system  mismatch,  we  used  all  SPAN  and  SCS930  positions  and  solved  for  the  affine 
transformation  (combination  of  translation,  rotation  and  dilation)  that  produced  the  best 
match,  in  a  least-squared  sense,  between  the  two  coordinate  systems.  The  rotation 
occurred  around  the  CP-2  control-point.  The  best-fitting  transformation  was  clockwise 
rotation  of  0.0440°,  corrected  scale  factors  of  1.001202  and  translation  of  3.0  cm  East,  9.6 
cm  North  and  -1769.5  cm  upwards.  The  affine  transformation  appears  to  have  removed 
any  bias  between  the  SPAN  and  SCS930  coordinate  systems  (Figure  23b). 

After  applying  the  affine  transformation,  we  then  compared  the  SPAN  predicted  position 
of  the  MT900  prism  against  the  position  measured  by  the  SCS930  for  the  straight  line 
and  S-curve  tests  (Figure  24)  and  the  left-  and  right-ramp  tests  (Figure  25).  RMS  errors 
and  standard  deviations  for  each  of  the  tests  are  given  in  Table  3.  For  the  straight-line 
path,  the  total  error  is  generally  less  than  1  cm,  with  RMS  errors  in  3-D  position  of 
between  0.78-0.93  cm.  Thus,  for  a  smooth,  straight-line  path,  the  Novatel  SPAN  is  able 
to  produce  errors  of  less  than  1  cm  in  position  and  better  than  0.4°  in  error  (recall  that 
with  the  lever-arm,  a  1  degree  error  in  orientation  results  in  a  2.5  cm  error  in  position). 
For  the  S-bends  there  is  a  larger  error  around  the  turns,  with  the  errors  increasing  to 
closer  to  2  cm.  Nevertheless,  the  3-D  RMS  errors  for  this  case  are  between  1.26  and  1.45 
cm.  When  the  left-hand  side  of  the  array  goes  over  a  ramp  there  is  a  relatively  large 
elevation  change  in  prism  position  and  a  corresponding  increase  in  the  positional  error 
(Figure  25a).  The  elevation  errors,  in  general,  appear  to  have  a  bias  in  this  run,  which 
may  indicate  a  systematic  error  in  the  estimated  GPS  elevation.  Three-dimensional  RMS 
errors  approach  3.25  cm,  but  the  horizontal  errors  are  smaller  with  RMS  less  than  1.8  cm. 
By  the  time  the  right-hand  side  ramp  test  is  conducted,  any  systematic  error  appears  to 
have  disappeared,  and  the  3-D  RMS  errors  are  between  1.5  and  2.5  cm. 


Figure  22:  Plan- view  (relative  to  CP-2)  of  the  prism  path  for  the  four  position  tests. 
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(a)  Direct  comparison  of  SCS  and  SPAN  positions 
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Figure  23:  Bias  in  the  SCS930  positions  before  (a)  and  after  (b)  an  affine  transformation. 
The  different  colors  represent  data  collection  events  with  the  towed-array  traveling  at 
different  speeds  (blue:  0.5;  red:  1.0;  green:  1.5;  cyan:  2.0;  and  magenta  3.0  mph). 
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(a)  Straight  path 


Distance  along  line  (m) 


(b)  Curved  path 


Distance  along  line  (m) 


Figure  24:  Difference  between  the  SPAN  predicted  position  of  the  MT900  prism  and  the 
SCS930  measured  position  of  the  prism  when  driving:  (a)  straight;  and  (b)  along  an  S-curve. 
Colors  are  different  speeds:  blue:  0.5;  red:  1.0;  green:  1.5;  cyan:  2.0;  and  magenta  3.0  mph. 


34 


September  2008 


MM-1571:  Next  Generation  Data  Collection  Platform 


Sky  Research,  Inc. 


(a)  Traversing  over  a  ramp  on  left  side  of  array 


Distance  along  line  (m) 


z  0  5  10  15  20  25  30  35 


Distance  along  line  (m) 


Distance  along  line  (m) 


(a)  Traversing  over  a  ramp  on  right  side  of  the  array 


Distance  along  line  (m) 


Distance  along  line  (m) 


Figure  25:  Difference  between  the  SPAN  predicted  position  of  the  MT900  prism  and  the 
SCS930  measured  position  of  the  prism  when  driving  over:  (a)  left-ramp;  and  (b)  right-ramp. 
Colors  are  different  speeds:  blue:  0.5;  red:  1.0;  green:  1.5;  cyan:  2.0;  and  magenta  3.0  mph. 
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Statistic 

Root  mean  sguared  error  (cm) 

Standard  deviation  (cm) 

Coordinate 

Easting 

Northing 

Elevation 

Total 

Easting 

Northing 

Elevation 

Total 

Straight  line 

0.57 

0.60 

0.49 

0.84 

0.47 

0.35 

0.31 

0.41 

0.22  m/s  (0.5  mph) 

0.43 

0.69 

0.65 

0.81 

0.37 

0.38 

0.38 

0.42 

0.9  m/s  (2  mph) 

0.64 

0.45 

0.25 

0.78 

0.45 

0.33 

0.25 

0.42 

1 .3  m/s  (3  mph) 

0.64 

0.67 

0.57 

0.93 

0.58 

0.32 

0.30 

0.39 

S-bends 

1.22 

0.58 

eeh 

0.50 

0.66 

0.56 

1.28 

0.65 

0.81 

mMmmssm 

1.11 

0.59 

KB 

0.46 

0.59 

0.67  m/s  (1.5  mph) 

1.36 

0.49 

0.57 

KE9 

0.47 

0.60 

0.9  m/s  (2  mph) 

1.16 

0.48 

0.52 

wm 

0.47 

0.50 

0.68 

1 .3  m/s  (3  mph) 

1.17 

0.79 

0.51 

mem 

1.14 

0.53 

0.41 

0.61 

Ramp  on  left 

1.01 

1.26 

2.17 

2.73 

0.43 

0.59 

0.62 

0.22  m/s  (0.5  mph) 

1.11 

1.35 

1.29 

2.18 

0.60 

0.53 

0.45  m/s  (1  mph) 

1.00 

1.31 

2.80 

3.25 

0.47 

0.41 

0.48 

0.48 

0.67  m/s  (1.5  mph) 

1.08 

1.45 

2.45 

3.04 

0.59 

0.43 

0.62 

0.70 

0.9  m/s  (2  mph) 

1.01 

1.09 

1.93 

2.44 

0.61 

0.49 

0.54 

0.69 

1 .3  m/s  (3  mph) 

0.84 

1.09 

2.37 

2.73 

0.63 

0.40 

0.72 

0.72 

Ramp  on  right 

0.66 

0.62 

1.64 

1.89 

0.62 

0.46 

0.52 

0.22  m/s  (0.5  mph) 

0.53 

0.58 

1.29 

1.51 

0.53 

0.45 

0.87 

0.62 

0.45  m/s  (1  mph) 

0.53 

0.46 

1.91 

2.03 

0.49 

0.31 

0.67 

0.68 

0.67  m/s  (1.5  mph) 

0.70 

0.72 

1.27 

1.62 

0.59 

0.44 

0.37 

0.36 

0.9  m/s  (2  mph) 

0.60 

0.67 

1.50 

1.75 

0.59 

0.53 

0.56 

1 .3  m/s  (3  mph) 

0.96 

0.69 

2.23 

2.52 

0.88 

0.40 

0.38 

Table  3:  Root  mean-squared  errors  and  standard  deviation  for  various  tests  on  the  smooth  surface. 


3.2.3  Rough  surface  tests 

We  repeated  the  smooth-road  tests  on  an  adjacent  field  where  the  towed  array  suffered 
larger  changes  in  orientation  and  elevation.  The  tests  were  conducted  on  the  second  day 
of  mobilization  and  utilized  a  new  set-up  of  the  SCS930  base-station.  An  affine 
transformation  was  again  required  to  remove  any  positional  biases.  Best  fitting 
parameters  were  clockwise  rotation  of  0.0853°,  corrected  scale  factors  of  1.00225  and 
translation  of  3.9  cm  East,  6.0  cm  North  and  -1771.0  cm  upwards.  The  results  of  the 
rough-surface  tests  (Figures  26  and  27,  Table  4)  are  similar  to  those  of  the  smooth  surface 
tests.  Several  of  the  datasets  display  biases  on  the  1  to  2  cm  level  on  the  reported 
elevations.  RMS  errors  on  the  straight-line  were  0.8  cm  on  the  smooth  surface  and  1.5  cm 
on  the  rough  surface,  for  the  S-curves  RMS  errors  were  1 .4  cm  on  both  surfaces,  for  the 
left  ramp  the  RMS  errors  were  2.7  cm  and  2.3  cm  and  for  the  right  ramp  they  were  1.9 
cm  and  1.5  cm  respectively.  Much  of  the  differences  in  the  smooth  and  rough-road 
scenarios  were  in  the  different  biases  in  the  elevation  angles.  Amongst  all  the  tests 
(smooth  and  rough  road,  different  speed,  straight,  s-bend,  ramps  etc)  the  RMS  horizontal 
error  varied  from  a  minimum  of  0.6  cm  to  a  maximum  of  1.7  cm.  These  represent  upper 
bounds  on  the  horizontal  positional  error  of  the  SPAN  system.  If  the  error  was  entirely 
due  to  an  error  in  the  SPAN  heading,  the  RMS  heading  error  would  vary  from  a 
minimum  of  0.2  0  to  a  maximum  of  0.6°.  These  results  are  remarkable  given  that  the 
SCS930  measured  positions  are  on  the  order  of  0.5  cm. 
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(a)  Straight  path 


Distance  along  line  (m) 


1-  0  5  10  15  20  25  30  35 


Distance  along  line  (m) 


(b)  Curved  path 


Figure  26:  For  the  rough  field  tests,  the  difference  between  the  SPAN  predicted  position  of 
the  MT900  prism  and  the  SCS930  measured  position  of  the  prism  when  driving:  (a)  straight; 
and  (b)  along  an  S-curve.  Colors  are  different  speeds:  blue:  0.5;  red:  1.0;  green:  1.5;  cyan: 
2.0;  and  magenta  3.0  mph 
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(a)  Traversing  over  a  ramp  on  the  left  side  of  the  array 


Distance  along  line  (m) 


z  0  5  10  15  20  25  30  35 


(b)  Traversing  over  a  ramp  on  the  right  side  of  the  array 


Figure  27:  For  the  rough  field  tests,  the  difference  between  the  SPAN  predicted  position  of 
the  MT900  prism  and  the  SCS930  measured  position  of  the  prism  when  driving  over:  (a)  left- 
ramp;  and  (b)  right  ramp.  Colors  are  different  speeds:  blue:  0.5;  red:  1.0;  green:  1.5;  cyan: 

2.0;  and  magenta  3.0  mph 
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Statistic 

Root  mean  squared  error  (cm) 

Standard  deviation  (cm) 

Coordinate 

Easting 

Northing 

Elevation 

Total 

Easting 

Northing 

Elevation 

Total 

Straight  line 

0.61 

0.82 

1.49 

1.88 

0.28 

0.61 

0.48 

0.22  m/s  (0.5  mph) 

0.45 

0.67 

0.89 

1.20 

0.20 

0.49 

0.36 

0.45  m/s  (1  mph) 

0.86 

1.12 

1.94 

2.40 

0.31 

0.32 

0.46 

0.52 

0.67  m/s  (1.5  mph) 

0.60 

0.87 

1.53 

1.86 

0.36 

0.30 

0.73 

0.49 

0.9  m/s  (2  mph) 

0.60 

1.13 

0.72 

1.47 

0.44 

0.31 

0.72 

0.39 

1 .3  m/s  (3  mph) 

0.54 

0.31 

2.38 

2.46 

0.54 

0.28 

0.65 

0.61 

S-bends 

0.90 

0.95 

1.37 

1.92 

0.73 

0.71 

0.61 

0.61 

0.22  m/s  (0.5  mph) 

0.67 

1.03 

1.34 

1.81 

0.66 

0.89 

0.51 

0.66 

0.45  m/s  (1  mph) 

0.78 

1.14 

1.49 

2.03 

0.75 

0.68 

0.48 

0.62 

0.67  m/s  (1.5  mph) 

0.81 

1.00 

0.99 

1.63 

0.69 

0.63 

0.77 

0.64 

0.9  m/s  (2  mph) 

0.98 

0.79 

1.16 

1.71 

0.72 

0.63 

0.61 

0.44 

1 .3  m/s  (3  mph) 

1.27 

0.81 

1.87 

2.41 

0.83 

0.71 

0.68 

0.66 

Ramp  on  left 

0.92 

0.37 

2.26 

2.54 

0.36 

0.69 

0.55 

0.22  m/s  (0.5  mph) 

0.60 

0.46 

2.97 

3.07 

0.48 

0.46 

0.49 

0.48 

0.45  m/s  (1  mph) 

0.63 

0.37 

2.97 

3.06 

0.42 

0.33 

0.70 

0.65 

0.67  m/s  (1.5  mph) 

0.69 

0.27 

2.18 

2.31 

0.45 

0.27 

0.84 

0.72 

0.9  m/s  (2  mph) 

1.27 

0.33 

1.99 

2.38 

0.48 

0.31 

0.59 

0.43 

1 .3  m/s  (3  mph) 

1.42 

0.42 

1.20 

1.90 

0.67 

0.41 

0.82 

0.50 

Ramp  on  right 

0.63 

0.40 

1.34 

1.56 

0.50 

0.30 

0.66 

0.50 

0.22  m/s  (0.5  mph) 

0.79 

0.64 

1.12 

1.51 

0.43 

0.26 

1.00 

0.54 

0.45  m/s  (1  mph) 

0.50 

0.42 

1.76 

1.87 

0.47 

0.41 

0.49 

0.49 

0.67  m/s  (1.5  mph) 

0.52 

0.31 

1.50 

1.61 

0.51 

0.30 

0.64 

0.52 

0.9  m/s  (2  mph) 

0.66 

0.32 

1.64 

1.80 

0.47 

0.30 

0.50 

1 .3  m/s  (3  mph) 

0.68 

0.33 

0.69 

1.02 

0.62 

0.25 

0.65 

Table  4:  Root  mean-squared  errors  and  standard  deviation  for  various  tests  on  the  rough  surface. 
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4  SIMULATIONS  OF  TEM  SYSTEMS 

In  this  section  of  the  report  we  describe  numerical  evaluation  of  different  candidate  sensor 
systems  designed  for  one-pass  detection  and  discrimination.  We  conducted  the  assessment 
using  both  semi-analytic  methods  (e.g.,  system  spatial  resolution,  model  uncertainty)  and 
monte-carlo  based  simulations.  The  semi-analytic  techniques  provide  powerful  metrics  to 
assess  the  expected  performance  of  different  systems,  but  do  not  easily  allow  positional  un¬ 
certainties  to  be  incorporated.  On  the  other-hand,  the  monte-carlo  simulations  can  incorpo¬ 
rate  all  types  of  errors  (including  positional  and  orientation  uncertainty)  but  the  results  can 
be  difficult  to  interpret  and  don’t  always  support  unambiguous  conclusions. 


4.1  CANDIDATE  SENSOR  SYSTEMS  AND  ASSESSMENT  PLAN 

All  of  the  candidate  sensor  systems  are  time-domain  system  with  equivalent  measurement 
characteristics  to  the  TEMTADS.  This  system  records  the  time-decay  in  a  user-configuration 
number  of  time-channels  between  40  /i  s  and  25  ms  after  pulse  turn-off.  ESTCP  Project 
MM0601  is  using  that  system  as  the  basis  for  a  cued-interrogation  mode  sensor  that  would 
collect  data  while  stationary  over  a  suspected  ordnance  item.  We  assume  the  time-gates  have 
been  configured  to  be  equivalent  to  the  EM-63  with  26  logarithmically  spaced  time-gates  be¬ 
tween  180  fis  and  25  ms  after  pulse  turn-off.  We  consider  three  different  array  configurations 
and  two  cart-based  configurations.  The  array-based  systems  are: 

1)  Cued-interrogation  array:  The  TEMTADS  array  with  25  transmitters  and  receivers  de¬ 
ployed  in  a  static  mode  (this  will  define  the  base-line  performance).  The  array  is  cen¬ 
tered  over  an  estimated  target  location,  and  then  transmitter  1  is  fired  several  times  with 
all  receivers  measuring  the  secondary  field.  Next  transmitter  2  is  fired  with  all  receivers 
measuring  the  data,  and  so-forth  until  all  25  transmitters  have  been  fired.  Figure  28  (a) 
is  a  schematic  view  of  25  transmitters  (red)/receivers  (blue),  where  each  loop  has  the 
size  of  0.4  m  x  0.4  m  and  the  receiver  loops  are  plotted  slightly  smaller  for  clear  display. 

2)  Modification  of  MTADS  cued-interrogation  array  to  conduct  one-pass  detection  and 
discrimination  by  replacing  the  25  transmitters  with  a  single  2  m  by  2  m  transmitter 
(but  keeping  the  25  receivers).  This  modified  system  is  shown  in  Figure  28  (b).  We 
will  refer  to  this  system  as  “Single-Tx/25-Rx”; 

3)  As  per  2)  but  with  four  1  m  by  1  m  transmitters  that  are  fired  sequentially  as  the  array 
transverses  over  an  area,  as  shown  in  28  (c).  We  will  refer  to  this  system  as  “Four- 
Tx/25-Rx”. 

These  two  modifications  to  the  cued  MTADS  array  represent  feasible  changes  to  the  system 
that  could  be  easily  implemented. 

We  consider  two  cart  systems: 

1)  Geonics  EM-63  equivalent  system:  with  a  0.5  m  by  0.5  m  receiver  coil  and  a  1  m  by  1 
m  transmitter. 
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Jt.  —  —  Jb  —  _  jl  —  —  -1 

(c)  (d) 


Figure  28:  Schematic  views  of  arrays  in  candidate  sensor  systems.  The  loop  in  red  represents 
a  transmitter,  in  blue  a  receiver,  (a)  The  cued-interrogation  array  of  25  transmitters  /receivers. 
Each  loop  is  0.4  m  x0.4  m.  (b)  The  array  of  single-Tx  (2  m  x2  m)/25-Rx.  (c)  Four-Tx/25- 
Rx,  transmitter:  1  m  x  1  m;  receiver:  0.4  m  x0.4  m.  (d)  EM63-M:  transmitter  of  1  m  x  1 
m.  Multi-component  receivers:  0.5  m  x  0.5  m. 


2)  Modification  of  the  EM-63  equivalent  to  include  horizontal  component  receivers  as 
shown  in  28  (d).  We  refer  to  this  modified  system  as  “EM63-M”,  where  M  indicates  a 
multi-component  receiver. 

The  EM-63  equivalent  will  provide  a  performance  base-line  of  the  ’’last  generation”  of  EMI 
sensors.  The  horizontal  coil  receiver  configurations  are  representative  of  next- generation 
systems  that  could  feasibly  be  deployed  in  a  one-pass  detection  and  discrimination  mode. 


4.2  SEMI-ANALYTICAL  APPROACH 


A  transient  EM  (TEM)  system  assumed  in  this  study  is  distributed  in  a  plane  above  the  surface 
to  interrogate  a  buried  object.  When  the  target  dimension  is  smaller  than  the  target-senor 
distance,  it  is  known  that  the  low  frequency  EMI  scattering  from  the  metallic  target  can  be 
well  described  by  an  equivalent  induced  dipole  (Bell.,  2005;  Das  et  al.,  1990;  Smith  et  al., 
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2004)  characterized  by  a  3  x  3  magnetic  polarizability  tensor 


Pit) 


Pn(t)  Pi2{t)  Piait) 
P2l(t)  P22(t)  P2ait) 
Pai{t)  P32{t)  Pzait) 


(6) 


The  elements  of  the  tensor  Pij(t)  represent  a  dipole  component  in  the  /th  Cartesian  direction 
due  to  a  primary  field  in  the  jth  Cartesian  direction  and  ptJ  =  pR.  The  equivalent  dipole  po¬ 
larizability  tensor  P(t)  contains  information  regarding  the  geometry  and  material  properties 
of  a  target  as  well  as  its  orientation. 

For  the  n-th  measurement  in  a  survey,  the  secondary  response  dn  at  time  instant  t  to  an 
equivalent  dipole  source  P(t)  located  at  r  can  be  written  as  (Das  et  al.,  1990) 

dnirRXn,t)  =  H^(r,  r&jP(<)Hr(r,  rTxJ,  (7) 

where  HR(r,  r^J  =  [HXR  HyR  HZ^\T  and  HT(r,rTxJ  =  [7/f  TR  Hf]T  are  the  magnetic 
fields  at  the  object  location  r  generated  by  the  receiver  and  transmitter  coils  with  their  centers 
at  rRxn  and  rTxn,  and  the  superscripts  x,  y,  and  z  denote  the  Cartesian  components  of  a  field 
and  the  superscript  T  denotes  a  transpose.  The  magnetic  fields  are  computed  via  line  integrals 
along  each  of  the  coils  (Das  et  al.,  1990).  Equation  (7)  describes  the  basic  EMI  process  of 
illuminating,  scattering  and  sensing  and  can  be  rearranged  in  another  scalar-product  form  as 

dnirRxn,t)  =  alir,rRxn,rTxn)q  (f),  (8) 

where  an(r,  r  Rxn .  rTxn)  is  a  6  x  1  column  vector  representing  spatial  sensitivities  of  the  n-th 
pair  of  transmitters-receivers  to  the  object  located  at  r,  and  q(f )  a  6  x  1  column  vector  whose 
components  are  the  elements  of  the  polarizability  tensor  P(t)  of  the  object: 


TJX  TJX 

tLRtLT 

'  Puit) " 

HxRHyT  +  HyRH*T 

Pn{t) 

TJX  TJZ  |  TJZ  TJX 

,q  it)  = 

Piait) 

HW 

P22(t) 

Km  +  KK 

P2a(t) 

TJZ  TJZ 

.  Paait)  _ 

With  a  total  of  N  measurements  the  data  and  model  parameter  relationship  can  be  conve¬ 
niently  expressed  in  vector-matrix  notation  as  follows, 

d(f)  =  A(r,rs)q(f)  +  n(f),  (10) 

where  d (t)  =  [di(t),  •  •  •  ,  d7v(t)]T  is  N  x  1  measured  data  vector  at  time  instant  t,  n (t)  the 
corresponding  N  x  1  additive  noise  vector,  and  /l(r.  rs)  is  a  N  x  6  matrix  denoting  the 
sensitivities  of  the  N  sensors  to  the  object  located  at  r, 


A(r,rs) 
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and  we  simply  condense  rs  =  [rRa,n,  rTxv\  to  represent  a  transmitter-receiver  configuration. 
Note  that  generally  N  =  Ns  x  NT  x  NR,  where  Ns  denotes  the  number  of  survey  locations, 
Nt  and  Nr  the  numbers  of  transmitters  and  receivers  designed  in  a  sensor  system.  Equation 
(10)  is  the  underlying  physical  model  used  for  semi-analytical  assessment  of  a  sensor  system. 


4.2.1  Linear  least-squares  solution 

Consider  an  EMI  survey  that  collects  data  with  a  specific  sensor  system  above  the  surface  to 
interrogate  a  buried  object.  The  objective  of  the  survey  is  to  estimate  the  target  polarization 
and  location  from  the  measured  data  through  the  relationship  defined  in  (10).  For  the  semi- 
analytical  approach,  we  consider  the  response  of  a  single-time  channel  data.  Thus  at  time  t, 
the  unknown  model  vector,  m,  is  given  by 


m  = 


(ID 


where  r  is  a  location  vector  of  an  object  and  q  is  composed  of  6  polarizabilities  at  time- 
channel  t  as  given  in  Equation  (9).  To  simplify  notation  we  will  stop  referring  to  the  explicit 
dependance  on  the  time  t.  The  inverse  problem  is  nonlinear  because  of  the  unknown  object 
location  and  its  coupling  with  the  unknown  dipolar  polarization,  although  the  latter  is  linearly 
related  to  the  measured  data,  as  observed  in  Equation  (10). 

For  a  nonlinear  problem,  the  observations  dQbs  can  be  generally  expressed  as 

dobs  =  ^[m]  +  n,  (12) 

where  T  denotes  a  forward  functional  operator  and  jE[m]  represents  the  theoretical  data 
which  are  computed  using  the  forward  model.  The  inverse  problem  can  be  stated  as  finding 
a  model  m  that  can  best  explain  the  observations  in  a  least-squares  sense,  i.e.,  minimize  the 
following  weighted  misfit, 

(j>d=\\Wd(dohs-F{m})\\2,  (13) 

where  W(i  is  a  weighting  matrix  (that  is  usually  related  to  the  inverse  covariance  matrix  of  the 
data). 

To  make  the  nonlinear  problem  tractable,  we  adopt  the  standard  approach  and  linearize  the 
problem  by  expanding  the  misfit  function  as  a  Taylors  series  in  the  vicinity  of  some  model 
m0.  We  retain  only  first-order  terms  so  that 


jE[m]  as  JE[m0]  +  J Am,  (14) 

where  Am  =  m  —  m0  is  a  model  correction  vector  and  J  is  N  x  M  the  Jacobian  matrix 
comprising  the  partial  derivatives  of  the  data  functional  evaluated  at  the  current  model, 


J 


dT 

dm 


m0> 


(15) 
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where  M  is  the  number  of  unknowns  which  is  9  in  our  case.  By  inserting  equation  (14)  into 
equation  (13),  we  obtain  a  linearized  approximation  to  the  misfit  function, 

0d=  ||Wd(Ad-  JAm)||2,  (16) 

where  Ad  =  dobs— J[m0]  are  the  data  residuals  between  the  measured  data  and  the  predicted 
ones  at  the  current  model  m0.  Equations.  (14)  through  (16)  define  a  linearized  problem. 
By  minimizing  Equation  (16),  we  obtain  the  weighted  least-squares  solution  to  the  model 
correction  vector  (Menke,  1989;  Tarantola,  1987)  at  one  step 

Am=(JTWjWdJ)~1JTWjWdAd.  (17) 

In  most  cases,  the  matrix  WjWd  in  Equation  (17)  can  be  chosen  as  the  inverse  of  the  covari¬ 
ance  matrix  Cd  of  data  errors  (if  it  can  be  estimated).  Using  the  notation  of  WjW,i  =  Cd 1 , 
Equation  (17)  can  be  simplified  slightly  as 

Am  =  (JTCJ1J)-1JTCJ1A  d.  (18) 

This  is  a  general  formal  solution  to  a  linear  or  linearized  problem.  If  errors  are  independent, 
Cd  is  a  diagonal  matrix  and  if  the  variance  of  each  observation  is  also  the  same,  then  each 
diagonal  element  of  Cd  will  be  equal  to  the  variance.  For  a  nonlinear  inverse  problem,  the 
linearization  procedure  in  Equations  (17)  or  (18)  may  be  applied  iteratively  to  attain  the 
optimal  estimate  of  the  model  parameters. 


4.2.2  Solution  uncertainty 

Uncertainty  in  the  estimated  model  parameters  can  be  assessed  in  terms  of  the  model  co- 
variance  matrix  (Menke,  1989;  Tarantola,  1987).  Denote  J+  as  the  generalized  inverse  in 
Equation  (18),  i.e.,  J+  =  (JT Cdl  J)-1  JT Cdl .  The  model  covariance  matrix  Cm  can  be 
written  as  (Menke,  1989;  Tarantola,  1987) 

Cm  =  J+CdJ+T  =  (JTC71  J)-1  (19) 

and  assuming  each  measurement  is  contaminated  by  Gaussian  white-noise  with  variance  cr^, 
then  Cd  =  crj/,  and 

Cm  =  cr2d(JTJ)-\  (20) 

The  model  parameter  covariance  matrix  depends  upon  the  uncertainty  in  the  data  as  well  as 
the  data  kernel  J  which  contains  the  sensitivities  of  the  data  functional  with  respect  to  the 
source  location  and  polarizability.  Thus  Cm  can  be  used  to  quantify  the  capability  of  a  sensor 
system  to  accurately  determine  the  polarization  tensor  and  location  parameter.  The  diagonal 
elements  of  Cm  contain  the  variance  of  each  model  parameter;  while  the  non-diagonal  ele¬ 
ments  describes  the  correlation  between  parameters  if  the  matrix  is  normalized.  The  sum  of 
the  diagonal  terms  of  the  covariance  matrix  can  be  used  as  a  single  measure  of  the  quality 
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of  a  given  survey.  For  example,  Liao  and  Carin  (2004)  used  the  trace  of  the  inverse  of  the 
Fisher  Information  Matrix  (FIM)  as  a  measure  of  survey  quality  in  the  development  of  an 
adaptive  framework  for  optimal  EMI  surveying  for  UXO  characterization.  The  inverse  of 
the  FIM,  is  known  as  the  Hessian  in  the  signal  processing  community  and  is  equivalent  to 
the  model  covariance  matrix.  In  the  linearized  problem  considered  here,  JTC^lJ  in  (19)  is  a 
Jacobian-based  approximation  to  the  Hessian,  which  is  often  used  when  solving  least-squares 
problems  (e.g.  Menke,  1989). 

Since  we  have  two  types  of  unknown  parameters  (the  object  location  and  its  polarizability), 
we  derive  two  measures  of  uncertainty  as  follows.  Partition  the  covariance  matrix  of  Cm  into 

Cm  =  (  •  •  •  Cp  )  ’  (21) 

where  Cr  and  Cp  are  block  matrices  of  size  3x3  and  6x6,  i.e.,  the  covariance  matrices 
corresponding  to  the  parameters  of  the  source  location  and  polarizability.  We  could  then,  for 
instance,  define  the  total  squared  uncertainty  of  the  polarizabilitties  as 

9 

=  Tr(C^  (22) 
4 

This  is  the  metric  that  was  used  by  Smith  et  al.  (2004)  to  assess  the  relative  merits  of  different 
transmitter/receiver  combinations  in  resolving  equivalent  dipole  polarizabilities. 

To  accurately  estimate  the  polarizability  parameters,  requires  a  well  determined  estimate  of 
the  source  location  r0.  For  the  two  sub-covariance  matrices  of  Cr  and  Cp,  we  can  derive 
analytical  expressions  to  determine  the  connection  between  the  different  parameters.  Without 
loss  of  generality,  set  Cd  =  /,  the  identity  matrix  (we  can  readmit  the  weighting  matrix  at  the 
end).  Divide  the  Jacobian  matrix  into  J  =  [D  A],  where  sub-matrices  D  and  A  express  the 
sensitivities  to  the  location  and  polarizabilty  of  an  object,  respectively.  We  can  then  expand 
the  LS  normal  equations  to  the  solution  (18)  as 

Dr  Ad  =  DtD  Ar  +  DT  AAq 
AT  Ad  =  ATD  Ar  +  ATAAq 

Consequently,  the  linearized  solution  of  Equation  (23)  can  be  written  as 

Ar  =  (DtD)~1Dt[ Ad  -  AAq] 

Aq  =  (GA^GAd  .  (24) 

G  =  At  -  AtD(DtD)~1Dt 

This  equation  is  equivalent  to  Equation  (18),  except  it  explicitly  separates  the  dependence  on 
location  from  that  on  polarizability.  Finally  the  model  uncertainties  in  Equation  (21)  for  both 
types  of  parameters  may  be  expressed  as 
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=  1(dtd)-1dt(i  -  a(ga)-'g\[(dtd)-1dt(i  -  Amr'er 

Cp  =  [{GA)-lG][{GA)-1G]T 

We  see  that  the  object  location  affects  the  estimation  of  polarizability  through  the  G  matrix, 
which  is  a  modification  of  A  by  incorporating  an  additional  term  that  involves  the  sensitivity 
matrix  D  to  the  object  location.  The  above  expressions  of  Equations  (24)  and  (25)  require  the 
inverse  of  the  matrices  DTD  and  GA.  The  conditioning  of  these  matrices  determines  how 
errors  in  the  data  translate  to  errors  in  the  location  and  polarizability. 


4.2.3  Survey-design  framework 

To  determine  the  effectiveness  of  a  given  sensor  system,  we  assume  a  survey  area  of  4  m 
x  4  m  with  an  object  located  at  the  area  center  and  buried  at  some  depth  (Figure  29).  The 
red  square  in  Figure  29  (a)  represents  the  cued  array  which  is  positioned  at  the  center  of  the 
region  and  occupies  an  area  of  2  m  x  2  m.  The  circles  represent  the  designated  stations  for 
the  moving  systems.  Figure  29  (b)  illustrates  how  a  4-Tx  system  is  sequentially  fired  when  it 
is  towed,  where  a  square  with  a  solid  line  represents  an  active  transmitter,  and  a  dashed  line 
represents  an  inactive  transmitter. 


Figure  29:  Schematic  view  of  the  experimental  surveys,  (a)  The  red  square  shows  the  area 
occupied  by  the  cued  array.  The  open  circles  represent  the  stations  for  the  towed  or  cart 
systems.  The  cross  represents  the  location  of  a  buried  object,  (b)  Schematic  view  of  4-Tx 
systems  fired  sequentially  when  moving.  A  square  with  a  solid  line  represents  an  active 
transmitter,  and  a  dashed  line  represents  an  inactive  transmitter. 


A  line-spacing  of  lx  =  1  m  is  used  for  the  multi-sensor  towed  system  and  lx  =  0.4  m  for 
the  cart  system  (to  provide  a  fair  comparison  as  the  receivers  in  the  towed  system  are  40  cm 
apart).  For  all  the  moving  systems,  the  station  spacing  is  sy  =  0.125  m  which  corresponds  to 
a  data  collection  rate  of  8  measurements  per  second  at  a  forward  speed  of  1  m/s.  The  sensors 
are  assumed  to  be  20  cm  above  the  surface. 

We  know  proceed  to  estimate  the  covariance  matrix  of  the  polarizabilities  and  the  location 
parameters  using  Equation  25.  We  first  start  with  a  spherical  test-object  and  then  consider  a 
test-object  with  cylindrical  symmetry  (representative  of  a  UXO). 
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4.2.4  Model  covariance  matrix  for  a  spherical  object 

We  first  assume  that  a  spherical  object  is  buried  at  (0,  0,  —0.4).  Because  it  is  a  sphere  it 
has  three  equal  principal  polarizations  and  we  assume  these  are  all  equal  to  9.81  (the  value 
of  the  principal  polarization  at  180  /is  for  a  37  mm  projectile).  Then  the  elements  of  the 
polarizability  tensor  P(t)  in  Equation  (6)  are  ptl  =  9.81  and  ptJ  =  0,1^  j.  Similar  to  Smith 
et  al.  (2004),  we  evaluate  the  Jacobian  J  at  the  true  parameter  values  so  as  to  remove  any 
dependence  of  algorithms  and  noise  realization  in  assessing  the  different  systems.  In  these 

Tables  5  list  the  covariance  matrices  for  equivalent  dipole  polarizabilities  computed  for  the  5 
sensor  systems  with  Equation  (25)  under  the  assumption  of  unit  variance  of  Cd  =  /.  Recall 
that  the  diagonal  elements  of  the  matrix  represent  variances  of  polarizability  elements  of  pn, 
P12,  P13,  P'22i  P23,  P33,  and  the  off-diagonal  elements  the  correlation  between  model  parameters 
if  the  matrix  is  diagonally  normalized. 

The  top  part  of  Table  5  shows  that  the  multi-Tx/Multi-Rx  cued  array  has  diagonal  covariance 
matrix  elements  on  the  order  of  10~3  with  the  off-diagonal  elements  even  smaller.  The  total 
polarizability  variance  (i.e.,  the  summation  of  the  diagonal  elements)  is  ep  =  0.023  .  In  Table 
5,  we  see  that  the  single-Tx/25-Rx  produces  small  values  of  variance  and  covariance,  with 
ep  =  0.2.  The  4-Tx/25-Rx,  has  slightly  better  accuracy  with  smaller  variances  in  model 
variances  and  ep  =  0.07.  The  multicomponent  EM-63  cart  system  cart  system  achieves 
similar  accuracy  to  the  two  one-pass  systems  with  ep  =  0.12.  In  contrast,  the  standard  EM- 
63  system  yields  the  largest  variance  in  polarizabilities  with  ep  =  10.34.  It  also  displays 
strong  correlations  between  pu,  p2 2,  P33,  with  the  correlation  coefficients  equal  to  0.97. 

For  N  independent  observations  the  uncertainty  is  proportional  to  1  /TV.  To  see  if  the  small 
uncertainties  in  the  new  systems  are  mainly  due  to  the  increased  number  of  data,  we  plot  the 
total  squared  polarizability  uncertainty  against  N  in  Figure  30.  On  the  plot,  we  also  show 
how  a2  /N  behaves  for  different  values  of  a2.  If  the  reduction  in  variance  was  simply  due 
to  an  increase  in  the  number  of  observations  then  different  system  would  lie  along  the  same 
a2/N  curve.  The  multi-sensor  systems  result  in  a  much  larger  reduction  in  variance  than 
1  /N,  indicating  that  they  have  inherently  superior  information  on  the  characteristics  of  the 
buried  object. 


4.2.5  Object  with  cylindrical  symmetry 

We  now  examine  the  changes  in  the  model  covariances  matrices  for  the  five  different  systems 
when  interrogating  an  object  with  cylindrical  symmetry  (UXO-like).  We  assume  a  polariza¬ 
tion  value  of  9.81  along  the  symmetry  axis  and  4.36  along  the  transverse  axis  (corresponding 
to  a  37  mm  projectile).  The  cylindrical  object  is  oriented  at  polar  angles  of  6  =  90°,  45°, 
and  0°,  corresponding  to  horizontal,  tilted,  and  vertical  scenarios,  respectively.  Note  that 
the  Jacobian  D  is  dependent  on  the  values  of  the  elements  of  and  thus  has  a  dependency 
on  the  orientation  of  the  object.  In  Table  6,  we  list  model  variances  of  locations  and  dipole 
polarizabilities  for  all  test  cases. 
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Number  of  observations 


Figure  30:  Comparison  of  total  squared  polarizability  uncertainty  against  the  number  of  ob¬ 
servations  for  a  sphere  at  0.4  m  depth.  The  solid  black  lines  show  how  a2 /N  behaves  for 
different  values  of  a2. 

For  the  horizontal  case,  the  6  polarizabilities  are  given  by 

\P11,P12,P13,P22,P23,P33]  =  [9.81,0.0,0.0,4.36,0.0,4.36]. 

The  polarizability  pu  is  related  to  the  axial  polarization,  while  p2 2  and  p33  are  related  to 
the  transverse  polarization.  All  systems  have  small  uncertainties  in  location  but  the  standard 
EM-63  produces  the  largest  values  for  locational  uncertainty.  Examining  the  polarizability 
uncertainty  results  in  Table  6  reveals  that  the  next-generation  systems  are  more  accurate  than 
the  standard  EM63.  The  EM63  has  the  largest  uncertainty  in  axial  polarization  with  relatively 
small  uncertainties  in  transverse  polarization.  Total  squared  uncertainty  for  the  EM63  is  3.07, 
significantly  larger  than  those  of  the  other  systems. 

When  the  object  is  tilted  at  45°,  the  polarizabilities  are 

[Pii,Pi2,Pi3,P22jP23,P33\  =  [7.09,0.0,-2.72,4.36,0.0,7.09]. 

where  pn  and  p33  are  related  to  the  average  of  the  axial  and  transverse  polarization,  p13  is 
half  the  difference  between  the  transverse  and  axial  polarizations,  while  p22  is  the  transverse 
polarization.  When  compared  with  the  horizontal  case,  the  magntitude  of  the  location  un¬ 
certainties  in  the  EM-63  are  increased.  This  is  illustrated  in  the  relevant  part  of  Table  6, 
where  the  corresponding  polarizability  uncertainties  (for  example  pn,  p33)  becomes  larger. 
Conversely,  the  modified  systems  still  produce  well  resolved  results.  In  this  case,  the  total 
squared  uncertainty  of  the  EM-63  is  8.55,  while  the  total  squared  uncertainties  of  the  other 
systems  remain  at  a  similar  order  to  the  horizontal  situation. 

In  the  vertical  case,  the  6  polarizabilities  are  given  by 

\PU,P12,P13,P22,P23,P33]  =  [4.36,0.0,0.0,4.36,0.0,9.81]. 

where  pu  and  p22  are  related  to  the  transverse  polarization,  while  p33  is  related  to  the  axial 
polarization  since  the  principal  axes  are  coincident  with  the  user  coordinate  axes.  This  is 
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the  worst  case  scenario  for  the  EM-63  as  indicated  by  the  very  large  uncertainty  in  the  axial 
polarization  in  Table  6.  The  new  systems  still  provide  very  satisfying  results. 

Overall,  the  new  systems  appear  promising  in  their  ability  to  accurately  resolve  the  polariz¬ 
abilities  of  an  object  at  various  orientations  and  different  depths.  In  terms  of  total  squared 
uncertainty,  the  system  of  4-Tx/25-Rx  is  the  best  one-pass  detection  candidate. 


4.3  FULL  NON-LINEAR  INVERSIONS:  EFFECTS  OF  MEASURE¬ 
MENTS  ERRORS  ON  DIPOLE  ESTIMATION 


We  now  wish  to  test  the  performance  of  the  candidate  systems  in  recovering  dipole  polar¬ 
izabilities  when  errors  in  sensor  orientation  and  location  are  included.  We  use  the  polariza¬ 
tion  tensor  model  in  a  parametric  form  in  which  the  i— th  principal  polarization  is  given  by 

Li(t)  =  (Pasion  and  Oldenburg,  2001).  The  model  parameters  m  contain  the  lo¬ 

cation  of  an  object,  its  orientation,  and  k,  /3,  7  parameters  along  three  principal  polarizations. 
Thus  for  a  three-polarization  inversion  of  EMI  data,  the  total  number  of  the  unknowns  is  15. 

To  produce  synthetic  data,  we  choose  three  items:  37  mm  projectile,  81  mm  mortar,  and 
4.2  inch  mortar,  representing  small,  medium,  and  large  UXOs  (with  polarizabilities  obtained 
from  previous  test-stand  measurements).  A  total  of  99  items  consisting  of  33  of  each  size  of 
UXO  are  randomly  seeded  over  cells  with  each  cell  taking  up  an  area  of  5  x5m2.  The  survey 
was  set  up  with  sensor  heights  of  h  =  0.20  m,  line  spacing  lx  —  1  m,  and  station  spacing 
ly  =  0.125  m.  Based  on  the  analysis  in  the  semi-analytical  approach,  a  line  spacing  of  0.4  m 
is  used  for  the  EM-63  for  a  fair  comparison.  The  synthetic  data  are  corrupted  with  simulated 
measurement  errors. 

Given  measured  data,  the  objective  is  to  find  the  m  that  minimizes  the  following  cost  function 
rh  =  argmin$(m) 

M  N 

$(m)  =  XX  tn ij  \dij0bs(tj)  r^j,pred (t j j  m)  |  (26) 

*= 1  j= 1 

subject  to  rrij  >  m  <  mu, 


where  irq  and  rn„  are  model  lower  and  upper  bounds  and  s  are  data  weighting  coefficients 
which  are  usually  chosen  as: 


1 

lJ  f  +  v\di,obs(tj)\' 


(27) 


The  parameter  7  is  a  percentage  that  down-weights  the  influence  of  the  large  data  values; 
typically  we  chose  1%  >  77  <  5%.  The  parameter  /  (ensures  that  the  denominator  in  the 
above  equation  does  not  get  too  close  to  zero  when  dhohs(tj)  approaches  zero. 


For  the  cued  array,  the  simulations  were  carried  out  under  assumption  of  Gaussian  random 
sensor  noise,  randn,  which  an  amplitude  as  function  of  time  of  — ^  *  randn.  We  use  a  =  0.6 
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for  the  dynamic  systems,  which  produces  noise  of  similar  magnitude  to  what  we’ve  observed 
at  field  sites  with  the  EM63.  The  cued  array  can  achieve  higher  signal-to-noise  ratio  (SNR) 
than  a  moving  system  through  stacking  (and  by  elimination  of  motion  noise).  For  the  cued  ar¬ 
ray,  we  ran  simulations  with  both  a  =  0.6  and  0.06  for  comparison  purposes.  For  the  moving 
systems,  the  simulations  were  done  considering  all  sources  of  data  error  that  typically  arise 
including  inaccuracies  in  sensor  position  and  orientation  as  well  as  sensor  noise.  We  assume 
that  the  errors  in  sensor  orientation  and  position  are  Gaussian  with  a  standard  deviation  of 
<r orien  =  2°,  a  standard  deviation  of  crp0s  =  2  cm,  and  5  cm. 

We  include  all  plots  and  tables  of  the  simulations  in  an  Appendix  to  demonstrate  and  describe 
the  recovered  parameters  and  polarizations  and  their  statistics.  In  the  main  text,  we  present  a 
summary  of  results  for  the  5  systems  in  a  few  combined  figures  and  tables  for  easy  examina¬ 
tion  and  refer  readers  to  the  Appendix  for  additional  details.  The  descriptions  or  definitions 
of  the  plots  appear  in  the  figures  that  are  given  in  the  appendix.  Briefly,  the  subscripts  p,  s, 
and  t  represent  the  primary,  secondary,  and  tertiary  parameter  or  polarization  in  the  analysis. 

Tables  7-9  list  all  SD  and  RMSE  (defined  in  Equations  (A-l)  and  (A-2))  values  to  give  a  total 
measure  of  difference  between  the  recovered  parameters  of  (k,  ft,  7)  and  the  true  parameters 
for  each  of  the  three  UXO  items.  In  the  summary  tables,  different  noisy  simulations  when 
using  the  moving  systems  refer  to: 

•  Noisy  case  I:  sensor  noise  of  *  randn  and  sensor  errors  with  crorien  =  2°,  crpos 
cm. 

•  Noisy  case  II:  sensor  noise  of  *  randn  and  sensor  errors  with  rrorien  =  2°,  <rpos 
cm. 

•  Noisy  case  III:  only  sensor  positional  errors  with  standard  error  of  <jpos  =  1  cm. 

•  Noisy  case  IV:  only  sensor  positional  errors  with  standard  error  of  crpos  =  2  cm. 

•  Noisy  case  V:  only  sensor  positional  errors  with  standard  error  of  crpos  =  5  cm. 

Figure  31  shows  the  summarized  results  for  the  5  systems  in  the  order  of  the  cued  array, 
single-Tx/25-Rx,  4-Tx/25-Rx,  and  EM63-M,  and  EM63  from  the  top  to  the  bottom.  For  the 
moving  systems  the  simulations  assumed  sensor  position  and  orientation  errors  of  crorien  =  2°, 
(j pos  =  2  cm.  In  the  kp-ks  plots,  we  generally  observe  three  tight  clusters  around  groundtruth 
values  for  the  first  4  systems  while  the  EM63  has  a  larger  spread  particularly  for  the  81  and 
4.2”  mortars.  In  the  (3P- plots,  the  cued  array  has  very  tight  parameter  clusters.  The  single- 
Tx,  4-Tx,  and  EM63-M  systems  have  three  linearly  clustered  patterns  but  the  EM63  has  no 
significant  clustering  for  the  ATC  81  mm.  In  the  plot  Lp(t\)  versus  Lp(t2i) / Lp(ti),  there  is 
a  tight-clustering  of  all  three  items  in  the  Cued,  4-Tx,  EM63-M,  and  Single-Tx  systems.  The 
single-Tx  system,  has  a  larger  spread  in  the  relative  decay  rate  for  the  37  mm  indicating  that 
the  polarization  at  late  time  is  not  well  recovered  (also  refer  to  Figures  A-3  in  the  appendix). 
For  the  EM63,  there  is  very  poor  clustering  for  any  of  the  items.  The  right  most  panel  shows 
the  difference  between  the  secondary  and  tertiary  polarizations  at  t 1  and  f2 1  normalized  by 
the  respective  true  values.  These  plots  may  be  used  to  examine  how  well  the  cylindrical 
symmetry  of  an  object  is  predicted.  Generally,  we  anticipate  that  polarization  at  early  time  is 
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recovered  more  accurately  than  that  at  late  time  and  thus  the  spread  is  most  likely  to  occur 
along  the  horizontal  direction  in  that  plot.  In  this  simulation,  we  can  see  that  the  EM63  has  a 
very  large  spread  that  would  prevent  use  of  this  feaature  vector  for  discrimination.  The  other 
systems  have  a  relatively  smaller  spread  thus  providing  an  indication  that  the  buried  item  has 
cylindrical  symmetry. 

Figure  32  presents  the  results  when  the  level  of  sensor  position  errors  is  increased  to  a  stan¬ 
dard  deviation  of  apos  =  5  cm  while  keeping  the  same  amount  of  sensor  noise  and  sensor- 
orientation  error.  The  towed  systems  of  single-Tx  and  4-Tx  still  display  good  performance. 
The  multi-component  EM63  cart  system  yields  acceptable  results  but  the  EM63  appears  quite 
sensitive  to  the  increasing  level  of  sensor  position  errors  and  produces  degraded  results. 

Overall,  the  simulations  show  that  the  4-Tx  towed  system  generally  has  smaller  SD  and 
RMSE  values  (as  summarized  in  tables  7-9)  and  achieves  the  best  performance  among  the 
moving  systems. 

We  also  conducted  simulations  for  the  EM63  under  pure  random  errors  in  sensor  position  to 
examine  what  minimal  level  of  positional  precision  is  justified  in  the  existing  EMI  systems. 
Three  levels  of  sensor  positional  errors  are  considered  in  the  simulation  with  standard  devi¬ 
ations  of  crpos  =  1  cm,  2  cm,  and  5  cm.  As  before,  we  put  all  the  plots  and  tables  derived 
from  these  simulation  in  the  appendix  and  summarize  the  results  in  Figure  33,  where  the 
results  obtained  previously  from  the  4-Tx  towed  system  are  also  displayed  for  comparison. 
For  small  positional  errors  of  crp0s  =  1  cm  and  2  cm,  the  EM63  produces  good  results  in 
which  all  three  items  can  be  distinguished  and  a  small  spread  in  the  discrepancy  between 
the  secondary  and  tertiary  polarizations  is  observed.  When  positional  errors  are  increased  to 
crpos  =  5  cm,  the  EM63  loses  the  capability  to  give  informative  results.  Figure  33  illustrates 
that  the  towed  4-Tx  system  with  a  positional  error  (Tpos  =  5  cm  has  similar  performance  to 
the  EM63  with  a  standard  deviation  of  less  than  <rpos  =  2  cm.  Notice  that  the  results  from 
the  4-Tx  system  were  obtained  considering  sources  of  errors  icluding  sensor  position  and 
orientation  as  well  as  sensor  noise.  Thus  we  can  expect  that  the  4-Tx  system  can  achieve  the 
desired  performance  without  the  strict  condition  of  1  cm  positional  accuracy  required  in  the 
existing  systems  for  successful  discrimination. 


4.4  SUMMARY  AND  DISCUSSION 

In  EMI  sensing,  our  ultimate  goal  is  to  accurately  recover  or  invert  for  parameters  of  a  buried 
target  that  allow  some  level  of  discrimination  to  be  achieved  between  UXOs  and  clutter.  The 
accuracy  of  the  inversion  solution  is  dependent  upon  a  number  of  factors  such  as  the  specifics 
of  the  inversion  method,  the  assumed  signal  model,  and  the  data  acquisition  attributes.  The 
sensor  system  used  is  probably  the  most  fundamental  determinant  of  discrimination  capabil¬ 
ity.  The  work  presented  here  was  focused  on  the  evaluation  of  next-generation  sensor  systems 
using  a  polarization  tensor  model. 

Methods  to  forecast  discrimination  performance  are  essential  for  sound  design  of  potential 
sensor  systems.  There  are  various  approaches  to  carry  out  such  a  task.  We  undertook  an 
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assessment  using  both  deterministic  and  stochastic  methods.  In  the  deterministic  framework, 
we  employed  semi-analytical  metrics,  based  on  a  linearization  of  the  underlying  inversion 
algorithm.  The  covariance  matrix  of  the  model  parameters  was  the  diagnostic  used  to  evalu¬ 
ate  candidate  systems.  Using  the  framework  we  predicted  that  the  multisensor  systems  had 
significantly  improved  performance  relative  to  the  EM-63.  The  4-Tx/25-Rx  system  was  the 
best  performing  of  the  one-pass  detection  and  discrimination  systems. 

In  the  Monte  Carlo  simulations  we  conducted  full  non-linear  inversion  using  synthetic  data 
that  are  contaminated  with  sensor  noise  and  errors  in  position  and  orienation.  In  the  simu¬ 
lations  we  considered  37  mm  projectiles,  and  81  and  4.2  inch  mortars.  The  tests  show  that 
the  standard  EM-63  requires  high  precision  of  position,  at  least  within  2  cm  accuracy  to  get 
reliable  estimates  of  dipole  polarizability.  They  further  showed  that  the  multi-sensor  systems 
perform  very  well.  In  particular,  the  4-Tx/25-Rx  system  is  robust  to  positional  and  orientation 
errors  of  5  cm  and  2°  and  outperforms  noise-free  EM-63  data  with  1-2  cm  positional  error. 
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Table  5:  Model  covariance  matrices  Cp  for  all  five  systems  for  the  the  sphere  at  z0  =  -0.4  m. 


Cued  array, 

ep  =  0.023 

model 

Pn 

Pi  2 

Pl3 

P22 

P23 

P33 

Pn 

0.007051 

0.000000 

0.000000 

0.003936 

0.000000 

0.004080 

Pi  2 

0.000000 

0.001558 

-0.000000 

0.000000 

0.000000 

0.000000 

P13 

0.000000 

-0.000000 

0.001047 

0.000000 

0.000000 

0.000000 

P22 

0.003936 

0.000000 

0.000000 

0.007051 

0.000000 

0.004080 

P23 

0.000000 

0.000000 

0.000000 

0.000000 

0.001047 

0.000000 

P33 

0.004080 

0.000000 

0.000000 

0.004080 

0.000000 

0.004889 

Single-Tx/25-Rx,  ep  =  0.2 

Pn 

0.068619 

0.000000 

0.000000 

0.017538 

0.000000 

0.023519 

Pl2 

0.000000 

0.021494 

-0.000000 

0.000000 

0.000000 

-0.000000 

Pl3 

0.000000 

-0.000000 

0.014334 

0.000000 

0.000000 

0.000000 

P22 

0.017538 

0.000000 

0.000000 

0.050879 

0.000000 

0.024003 

P23 

0.000000 

0.000000 

0.000000 

0.000000 

0.011604 

0.000000 

P33 

0.023519 

-0.000000 

0.000000 

0.024003 

0.000000 

0.035774 

4-Tx/25-Rx,  ep  =  0.07 

Pn 

0.027855 

-0.000232 

-0.000054 

0.007994 

0.000048 

0.008712 

P\2 

-0.000232 

0.005428 

0.000041 

-0.000140 

0.000054 

-0.000168 

Pl3 

-0.000054 

0.000041 

0.003439 

-0.000022 

-0.000083 

-0.000068 

P22 

0.007994 

-0.000140 

-0.000022 

0.016335 

0.000085 

0.010810 

P23 

0.000048 

0.000054 

-0.000083 

0.000085 

0.002846 

0.000185 

P33 

0.008712 

-0.000168 

-0.000068 

0.010810 

0.000185 

0.012823 

EM63-M, 

ep  =  0.12 

Pn 

0.017193 

-0.000000 

-0.005574 

0.005462 

0.000000 

0.013728 

Pl2 

-0.000000 

0.003046 

-0.000000 

-0.000000 

-0.000564 

-0.000000 

Pl3 

-0.005574 

-0.000000 

0.024482 

0.002195 

0.000000 

0.002832 

P22 

0.005462 

-0.000000 

0.002195 

0.013430 

0.000000 

0.012876 

P23 

0.000000 

-0.000564 

0.000000 

0.000000 

0.018842 

0.000000 

P33 

0.013728 

-0.000000 

0.002832 

0.012876 

0.000000 

0.039962 

EM63,  ep 

=  10.34 

Pn 

2.021576 

0.000000 

0.096757 

1.992680 

-0.000000 

3.439111 

Pl2 

0.000000 

0.032325 

0.000000 

0.000000 

0.002251 

0.000000 

Pl3 

0.096757 

0.000000 

0.089334 

0.119573 

0.000000 

0.207360 

P22 

1.992680 

0.000000 

0.119573 

2.085530 

0.000000 

3.500249 

P23 

-0.000000 

0.002251 

0.000000 

0.000000 

0.074910 

0.000000 

P33 

3.439111 

0.000000 

0.207360 

3.500249 

0.000000 

6.038188 
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Table  6:  Model  variances  for  the  object  with  cylindrical  symmetry. 


System 

< 

*£  bO 

< 

<2 

<3 

2 

<JP22 

2 

°P23 

<33 

ep 

9  =  90°,  z0  = 

-0.4  m,  [pii,pi2,Pi3,P22,P23,P33]  =  [9.81,0.0,0.0,4.36,0.0,4.36]. 

Cued  array 

3.0E-6 

4.0E-6 

4.0E-6 

0.018 

1.6E-3 

7.3E-4 

6.0E-3 

8.8E-4 

3.8E-3 

0.031 

Single  Tx 

4.9E-5 

5.7E-5 

4.9E-5 

0.12 

0.02 

0.007 

0.044 

0.010 

0.026 

0.22 

4  Txs 

1.2E-5 

1.6E-5 

1.6E-5 

0.049 

0.0054 

0.0018 

0.015 

0.0024 

0.011 

0.085 

EM63-M 

1.6E-5 

1.5E-5 

1.4E-5 

0.043 

0.0031 

0.0035 

0.012 

0.007 

0.023 

0.091 

EM63 

9.0E-5 

9.9E-5 

7.4E-4 

1.90 

0.033 

0.025 

0.27 

0.039 

0.81 

3.17 

9  =  45u 

,  Z0  =  -0.4  m,  [pn,  Pi2,  Pl3,  P22,  P23,  P33} 

=  [7.09,  0.0,  -2.72, 4.36,  0.0,  7.09], 

Cued  array 

3.0E-6 

4.0E-6 

2.0E-6 

8.4E-3 

1.8E-3 

1.8E-3 

4.6E-3 

1.5E-3 

4.8E-3 

0.023 

Single  Tx 

3.7E-5 

3.7E-5 

2.6E-5 

0.074 

0.024 

0.019 

0.035 

0.015 

0.036 

0.203 

4  Txs 

9.0E-5 

1.0E-5 

7.0E-6 

0.029 

0.0061 

0.005 

0.009 

0.004 

0.013 

0.067 

EM63-M 

1.7E-5 

1.5E-5 

1.0E-5 

0.042 

0.0078 

0.018 

0.0094 

0.015 

0.038 

0.13 

EM63 

2.3E-4 

2.1E-4 

1.9E-3 

2.35 

0.074 

1.05 

0.32 

0.163 

4.60 

8.55 

9  =  0 

0  „  _ 

>  Zo  -  - 

0.4  m,  \pu 

)  Pl2i  Pl3  j  P22  j  P23  j  P33 

=  [4.36,0.0,0.0,4.36,0.0,9.81], 

Cued  Array 

4.0E-6 

4.0E-6 

1.0E-6 

4.0E-3 

1.6E-3 

2.7E-3 

4.0E-3 

2.7E-3 

5.6E-3 

0.021 

Single  Tx 

2.8E-5 

2.8E-5 

1.8E-5 

0.054 

0.021 

0.024 

0.033 

0.022 

0.044 

0.20 

4  Txs 

8.0E-6 

8.0E-6 

4.0E-6 

0.023 

5.4E-3 

6.8E-3 

7.8E-3 

6.5E-3 

0.014 

0.064 

EM63-M 

2.9E-5 

1.3E-5 

8.0E-6 

0.015 

3.0E-3 

0.048 

0.010 

0.024 

0.051 

0.15 

EM63 

5.0E-4 

4.8E-4 

2.20E-3 

0.22 

0.032 

0.68 

0.20 

0.66 

10.38 

12.16 

9  =  0 

u,  z0  =  - 

0.2  m.  [pu,Pl2,Pl3,P22,P23,P33 

=  [4.36,0.0,0.0,4.36,0.0,9.81], 

Cued  array 

0.00 

0.00 

0.00 

2.4E-4 

1.1E-4 

1.2E-4 

2.4E-4 

1.2E-4 

3.0E-4 

0.001 

Single  Tx 

0.00 

0.00 

0.00 

0.012 

2.5E-3 

3.2E-3 

3.5E-3 

1.8E-3 

4.0E-3 

0.027 

4  Txs 

0.00 

0.00 

0.00 

3.8E-3 

4.2E-4 

5.7E-4 

5.5E-4 

3.9E-4 

9.8E-4 

6.7E-3 

EM63-M 

0.0 

0.0 

0.0 

6.2E-4 

1.56E-4 

8.5E-4 

3.3E-4 

1.56E-4 

3.7E-4 

2.5E-3 

EM63 

2.0E-5 

1.6E-5 

1.3E-4 

0.025 

0.0023 

0.04 

0.017 

0.03 

0.76 

0.88 
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Figure  31:  The  results  of  three-polarization  inversion  of  noisy  data  using  5  systems.  Five 
rows:  (a)  Cued  array,  sensor  noise  only;  (b)  Single-Tx/25-Rx;  (c)  4-Tx/25-Rx;  (d)  EM63- 
M;  (e)  EM63.  For  (b)-(e)  measurement  errors:  sensor  noise  +  errors  in  sensor  position  and 
orientation  with  crorien  =  2°,  <jpos  =  2  cm.  In  all  cases,  sensor  noise  is  simulated  as  *randn. 
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Figure  32:  The  results  of  three-polarization  inversion  of  noisy  data  using  5  systems.  Five 
rows:  (a)  Cued  array,  sensor  noise  only;  (b)  Single-Tx/25-Rx;  (c)  4-Tx/25-Rx;  (d)  EM63- 
M;  (e)  EM63.  For  (b)-(e)  measurement  errors:  sensor  noise  +  errors  in  sensor  position  and 
orientation  with  crorien  =  2°,  <rpos  =  5  cm.  In  all  cases,  sensor  noise  is  simulated  as  *randn. 
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Figure  33:  The  results  of  three-polarization  inversion  of  noisy  data  using  EM63  and  4-Tx 
systems.  For  4-Tx  system,  measurement  errors  including  sensor  noise  (^|  *  randn),  pertur¬ 
bations  in  sensor  position  and  orientation  with  crorien  =  2°,  apos  =  2,  5  cm  in  row  (a)  and  (b). 
For  EM63:  only  perturbations  in  sensor  position  are  simulated  with  erpos  =  1,  2,  5  cm  in  row 
(c)-(e). 
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Table  7:  Statistical  results  for  3 -polarization  inversion  of  37  mm  using  5  systems. 


Parameters 

hp 

ks 

h 

/ % 

Ps 

Pt 

7p 

7s 

It 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

SD 

0.06 

0.01 

0.01 

0.00 

0.00 

0.01 

0.01 

0.33 

0.25 

Cued  array,  noise:  ^=7 

RMSE 

0.01 

0.01 

0.01 

0.00 

0.00 

0.00 

0.00 

0.06 

0.04 

SD 

0.26 

0.03 

0.06 

0.02 

0.02 

0.04 

0.28 

2.04 

3.26 

Cued  array,  noise:  . 

RMSE 

0.05 

0.04 

0.09 

0.03 

0.02 

0.03 

0.07 

0.36 

0.57 

Single-Tx,  noisy  case  I 

SD 

0.28 

0.09 

0.18 

0.03 

0.07 

0.13 

2.37 

6.20 

8.92 

RMSE 

0.06 

0.14 

0.29 

0.06 

0.06 

0.11 

0.56 

1.09 

1.57 

4-Tx,  noisy  case  I 

SD 

0.26 

0.06 

0.08 

0.01 

0.05 

0.06 

0.38 

5.08 

2.73 

RMSE 

0.05 

0.10 

0.12 

0.03 

0.04 

0.05 

0.09 

0.90 

0.48 

EM63-M,  noisy  case  I 

SD 

0.58 

0.12 

0.28 

0.04 

0.15 

0.16 

0.53 

4.62 

6.23 

RMSE 

0.12 

0.19 

0.44 

0.09 

0.13 

0.14 

0.12 

0.81 

1.10 

EM63,  noisy  case  I 

SD 

1.45 

0.16 

0.29 

0.05 

0.14 

0.18 

3.49 

4.73 

8.76 

RMSE 

0.30 

0.26 

0.46 

0.11 

0.12 

0.16 

0.82 

0.83 

1.54 

Single-Tx,  noisy  case  II 

SD 

0.51 

0.19 

0.31 

0.03 

0.10 

0.17 

2.52 

5.54 

10.20 

RMSE 

0.11 

0.30 

0.50 

0.07 

0.09 

0.15 

0.59 

0.98 

1.80 

4-Tx,  noisy  case  II 

SD 

0.73 

0.12 

0.14 

0.02 

0.08 

0.12 

0.64 

8.67 

9.87 

RMSE 

0.15 

0.19 

0.23 

0.05 

0.07 

0.11 

0.15 

1.53 

1.74 

EM63-M,  noisy  case  II 

SD 

1.79 

0.22 

0.41 

0.07 

0.22 

0.28 

1.89 

9.87 

11.35 

RMSE 

0.37 

0.36 

0.65 

0.17 

0.19 

0.24 

0.44 

1.74 

2.00 

EM63,  noisy  case  II 

SD 

2.64 

0.31 

0.36 

0.08 

0.25 

0.29 

2.69 

9.58 

7.62 

RMSE 

0.55 

0.50 

0.58 

0.17 

0.22 

0.25 

0.63 

1.69 

1.34 

EM63,  noisy  case  III 

SD 

0.59 

0.07 

0.08 

0.01 

0.03 

0.04 

0.09 

1.96 

2.78 

RMSE 

0.12 

0.11 

0.14 

0.02 

0.02 

0.04 

0.02 

0.35 

0.49 

EM63,  noisy  case  IV 

SD 

1.41 

0.16 

0.22 

0.02 

0.06 

0.10 

0.19 

4.28 

7.85 

RMSE 

0.29 

0.25 

0.35 

0.04 

0.05 

0.09 

0.04 

0.75 

1.38 

EM63,  noisy  case  V 

SD 

2.43 

0.33 

0.42 

0.06 

0.29 

0.26 

4.27 

6.33 

12.74 

RMSE 

0.50 

0.53 

0.67 

0.14 

0.25 

0.23 

1.00 

1.12 

2.25 
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Table  8:  Statistical  results  for  3 -polarization  inversion  of  81  mm  using  5  systems. 


Parameters 

kp 

ks 

h 

Pp 

Ps 

Pt 

7P 

7s 

7t 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

SD 

0.04 

0.05 

0.07 

0.00 

0.00 

0.00 

0.02 

0.02 

0.02 

Cued  array,  noise: 

RMSE 

0.00 

0.01 

0.01 

0.00 

0.00 

0.00 

0.01 

0.01 

0.01 

SD 

0.80 

0.59 

0.48 

0.01 

0.02 

0.03 

0.26 

0.22 

0.22 

Cued  array,  noise: 

RMSE 

0.03 

0.06 

0.05 

0.02 

0.03 

0.04 

0.07 

0.12 

0.12 

Single-Tx,  noisy  case  I 

SD 

1.95 

1.42 

1.27 

0.03 

0.06 

0.07 

1.60 

4.84 

4.84 

RMSE 

0.08 

0.15 

0.13 

0.05 

0.08 

0.09 

0.43 

2.59 

2.59 

4-Tx,  noisy  case  I 

SD 

1.54 

1.17 

0.83 

0.01 

0.06 

0.05 

0.65 

2.20 

2.20 

RMSE 

0.06 

0.12 

0.08 

0.02 

0.07 

0.06 

0.17 

1.18 

1.18 

EM63-M,  noisy  case  I 

SD 

4.45 

1.80 

2.04 

0.05 

0.09 

0.08 

8.38 

1.71 

1.71 

RMSE 

0.19 

0.19 

0.21 

0.08 

0.11 

0.10 

2.24 

0.92 

0.92 

EM63,  noisy  case  I 

SD 

21.48 

3.59 

13.15 

0.07 

0.20 

0.15 

8.11 

8.71 

8.71 

RMSE 

0.90 

0.37 

1.35 

0.11 

0.25 

0.20 

2.17 

4.65 

4.65 

SD 

2.71 

1.88 

1.69 

0.03 

0.07 

0.09 

2.83 

4.86 

4.86 

Single-Tx,  noisy  case  II 

RMSE 

0.11 

0.19 

0.17 

0.05 

0.08 

0.11 

0.76 

2.60 

2.60 

4-Tx,  noisy  case  II 

SD 

2.52 

1.62 

1.54 

0.03 

0.08 

0.08 

1.50 

1.74 

1.74 

RMSE 

0.11 

0.17 

0.16 

0.06 

0.10 

0.10 

0.40 

0.93 

0.93 

EM63-M,  noisy  case  II 

SD 

7.28 

2.42 

6.06 

0.04 

0.14 

0.12 

3.30 

6.22 

6.22 

RMSE 

0.31 

0.25 

0.62 

0.06 

0.18 

0.15 

0.88 

3.32 

3.32 

EM63,  noisy  case  II 

SD 

34.96 

4.69 

9.74 

0.07 

0.13 

0.13 

3.78 

4.39 

4.39 

RMSE 

1.46 

0.48 

1.00 

0.11 

0.17 

0.16 

1.01 

2.35 

2.35 

EM63,  noisy  case  III 

SD 

3.95 

0.99 

1.27 

0.01 

0.01 

0.01 

0.10 

0.22 

0.22 

RMSE 

0.17 

0.10 

0.13 

0.01 

0.02 

0.01 

0.03 

0.12 

0.12 

EM63,  noisy  case  IV 

SD 

7.22 

2.44 

2.28 

0.02 

0.03 

0.06 

3.13 

1.44 

1.44 

RMSE 

0.30 

0.25 

0.23 

0.03 

0.04 

0.07 

0.84 

0.77 

0.77 

EM63,  noisy  case  V 

SD 

31.56 

4.65 

5.27 

0.03 

0.14 

0.05 

0.56 

2.09 

2.09 

RMSE 

1.32 

0.48 

0.54 

0.05 

0.17 

0.06 

0.15 

1.12 

1.12 
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Table  9:  Statistical  results  for  3 -polarization  inversion  of  4.2  inch  mortars  using  5  systems. 


Parameters 

kp 

hs 

kt 

A j 

Pa 

Pt 

7 p 

7s 

7 1 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

SD 

0.22 

0.25 

0.14 

0.00 

0.01 

0.01 

0.02 

0.16 

0.07 

Cued  array,  noise:  ^=7 

RMSE 

0.00 

0.01 

0.01 

0.00 

0.01 

0.01 

0.01 

0.04 

0.02 

SD 

2.05 

0.68 

0.56 

0.01 

0.01 

0.02 

0.17 

0.27 

0.32 

Cued  array,  noise: 

RMSE 

0.02 

0.03 

0.02 

0.01 

0.01 

0.02 

0.04 

0.07 

0.08 

Single-Tx,  noisy  case  I 

SD 

8.53 

2.09 

1.61 

0.03 

0.03 

0.03 

0.55 

2.69 

0.83 

RMSE 

0.09 

0.09 

0.07 

0.05 

0.03 

0.04 

0.13 

0.71 

0.22 

4-Tx,  noisy  case  I 

SD 

5.41 

2.45 

1.20 

0.02 

0.04 

0.04 

0.43 

3.82 

0.56 

RMSE 

0.05 

0.11 

0.05 

0.04 

0.04 

0.04 

0.10 

1.01 

0.15 

EM63-M,  noisy  case  I 

SD 

8.30 

3.29 

3.48 

0.05 

0.07 

0.08 

1.76 

4.09 

0.88 

RMSE 

0.08 

0.15 

0.15 

0.10 

0.07 

0.09 

0.41 

1.08 

0.23 

EM63,  noisy  case  I 

SD 

31.76 

6.86 

12.32 

0.04 

0.09 

0.12 

1.35 

7.35 

5.38 

RMSE 

0.32 

0.30 

0.55 

0.07 

0.10 

0.13 

0.31 

1.93 

1.42 

Single-Tx,  noisy  case  II 

SD 

14.14 

3.72 

2.18 

0.03 

0.05 

0.06 

0.43 

2.42 

1.45 

RMSE 

0.14 

0.16 

0.10 

0.06 

0.05 

0.06 

0.10 

0.64 

0.38 

4-Tx,  noisy  case  II 

SD 

13.77 

4.06 

2.16 

0.05 

0.10 

0.07 

0.35 

3.35 

0.86 

RMSE 

0.14 

0.18 

0.10 

0.09 

0.11 

0.08 

0.08 

0.88 

0.23 

EM63-M,  noisy  case  II 

SD 

23.91 

6.20 

6.03 

0.06 

0.10 

0.12 

3.43 

8.49 

1.78 

RMSE 

0.24 

0.27 

0.27 

0.12 

0.10 

0.12 

0.79 

2.23 

0.47 

EM63,  noisy  case  II 

SD 

44.10 

8.03 

10.35 

0.06 

0.15 

0.12 

10.44 

11.44 

1.56 

RMSE 

0.45 

0.36 

0.46 

0.11 

0.15 

0.12 

2.41 

3.01 

0.41 

EM63,  noisy  case  III 

SD 

8.83 

2.16 

1.46 

0.01 

0.01 

0.01 

0.05 

0.24 

0.40 

RMSE 

0.09 

0.10 

0.06 

0.01 

0.01 

0.01 

0.01 

0.06 

0.11 

EM63,  noisy  case  IV 

SD 

23.94 

3.94 

7.85 

0.01 

0.03 

0.03 

0.12 

0.73 

0.83 

RMSE 

0.24 

0.17 

0.35 

0.02 

0.03 

0.04 

0.03 

0.19 

0.22 

EM63,  noisy  case  V 

SD 

142.04 

8.91 

23.03 

0.09 

0.14 

0.13 

0.80 

2.63 

5.27 

RMSE 

1.45 

0.39 

1.02 

0.17 

0.14 

0.13 

0.18 

0.69 

1.39 

60 


September  2008 


MM-1571:  Next  Generation  Data  Collection  Platform 


Sky  Research,  Inc. 


5  DISCUSSION  AND  RECOMMENDATIONS 

5.1  DISCUSSION 

In  this  project,  we  conducted  a  feasibility  analysis  of  a  one-pass  detection  and 
discrimination  system.  We  focused  on  (i)  estimating  the  position  and  orientation 
accuracies  that  can  be  achieved  in  a  dynamic  survey  mode;  and  (ii)  on  simulating  the 
discrimination  performance  of  a  number  of  candidate  EMI  sensors  when  subjected  to 
different  position  and  orientation  errors. 

For  sensor  positioning  we  used  Novatel  SPAN  system  with  a  tactical  grade  IMU  that  has 
a  drift  rate  of  1°  hr"1.  The  data  were  post-processed  in  the  commercially  available  Inertial 
Explorer  package.  Two  positional  tests  were  conducted,  the  first  using  a  Rotating 
Amusement  Park  carousel,  and  the  second  using  an  EM61  towed  array. 

The  RAP  was  instrumented  with  a  precision  optical  encoder  and  the  SPAN  system  was 
installed  on  a  tilt-table  at  a  fixed  position  on  the  RAP.  The  known  distance  of  the  SPAN 
system  from  the  axis  of  revolution  of  the  RAP  ride,  the  tilt  table  angle  and  the  optical 
encoder  information  provided  accurate  ground-truth.  The  GPS  and  IMU  data  from  the 
SPAN  were  augmented  by  an  optical  encoder  attached  to  a  bicycle  wheel  that  was 
attached  to  the  RAP  and  which  was  in  constant  contact  with  the  ground.  For  one  run, 
when  the  system  was  moving  at  a  brisk  2  meters  per  second,  the  position  estimated  from 
the  IMU/GPS  combination  had  a  root-mean-square  error  of  between  1.1  and  1.8  cm. 

For  the  second  position  test,  a  laser  prism  was  mounted  on  the  rear  left  comer  of  an 
EM61  towed-array,  about  1.5  m  away  from  the  Novatel  SPAN  system  which  was 
mounted  on  the  tow-bar  in  line  with  the  center  of  the  towed-array.  The  SPAN  predicted 
positions  of  the  prism  were  compared  to  the  prism  positions  measured  with  a  Trimble 
SCS930  Universal  Total  Station.  Data  were  collected  at  various  speeds  between  0.2  and 
1 .3  meters  per  second,  on  both  smooth  and  rough  ground,  with  the  array  traveling  along, 
straight  or  curved  paths  and  with  ramps  sometimes  present  that  one  side  of  the  cart  had  to 
travel  over.  Amongst  all  the  tests  the  root-mean-squared  (RMS)  horizontal  error  varied 
from  a  minimum  of  0.6  cm  to  a  maximum  of  1.7  cm.  These  represent  upper  bounds  on 
the  horizontal  positional  error  of  the  SPAN  system.  If  the  error  was  entirely  due  to  the 
SPAN  heading,  the  RMS  heading  error  would  vary  from  a  minimum  of  0.2°  to  a 
maximum  of  0.6°.  RMS  errors  of  the  3-D  positions  (including  elevation)  varied  from  0.8 
to  3.2  cm.  The  runs  with  larger  errors  typically  exhibited  systematic  biases  in  the 
elevation  angles.  This  indicates  that  there  may  be  1-2  cm  magnitude  systematic 
differences  in  the  SPAN  elevations  along  adjacent  transect  paths. 

Numerical  investigation  of  simulation  performance  occurred  along  two  fronts.  The  first 
was  semi-analytical  and  based  on  estimating  the  expected  errors  in  the  positions  and 
polarizabilities  of  a  buried  metallic  object.  The  second  was  based  on  Monte-Carlo 
simulations  of  systems  subjected  to  different  amounts  of  sensor,  position  and  orientation 
noise.  Five  different  sensor  systems  were  evaluated  including  an  EM63  equivalent,  EM- 
63  with  3-axis  receiver  and  three  versions  of  the  TEMTADS  with  one  deployed  in  cued 
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mode,  and  the  other  two  in  dynamic  mode.  Both  the  semi-analytical  and  monte-carlo 
approaches  demonstrated  the  excellent  performance  of  the  TEMTADS  cued-system. 
They  also  revealed  that  the  dynamic  variants  of  that  system  have  significantly  improved 
performance  over  an  EM63  equivalent.  In  particular,  the  multiple  receiver  coils  makes 
the  TEMTADS  more  tolerant  of  positional  errors  the  EM63.  The  four-transmitter 
TEMTADS  with  5  cm  position  error  and  2°  orientation  error  outperforms  an  EM63  with 
1  cm  and  1°  error. 

5.2  RECOMMENDATIONS 

We  recommend  adapting  the  TEMTADS  cued-interrogation  array  for  use  in  a  one-pass 
detection  and  discrimination  mode.  The  following  issues  need  to  be  considered  when 
adapting  the  TEMTADS  for  use  in  a  one-pass  detection  and  discrimination  mode: 

1)  Positional  and  orientation  accuracy:  In  a  dynamic  mode,  discrimination  relies  on 
accurate  relative  position  between  successive  measurements.  The  simulations  in  the 
last  section  demonstrated  that  parameter  clusters  produced  with  one  suggested  mobile 
version  of  the  TEMTADS  deployed  with  5  cm  error  in  position  and  2  degree  error  in 
orientation  were  comparable  to  clusters  produced  from  EM63  data  with  1  cm  position 
and  1  degree  orientation  error.  Furthermore,  the  tests  with  the  EM61  towed  array  at 
FLBGR  demonstrated  SPAN  RMS  position  errors  of  better  than  1-2  cm  and 
orientations  of  better  than  0.6  degrees.  We  therefore  recommend  using  the  Novatel 
SPAN  system  with  data  post-processed  in  the  Inertial  Explorer  software  package, 
with  the  GPS  potentially  augmented  with  the  Trimble  SCS930  laser-positioning 
system. 

2)  Timing  accuracy:  The  accurate  and  stable  synchronization  of  the  different  sensor 
components  is  a  critical  requirement  for  successful  discrimination  based  on 
dynamically  collected  data.  We  recommend  using  the  SKY  DAS  to  provide  accurate 
time-stamping  of  all  system  components.  It  has  the  ability  to  simultaneously  time- 
stamp  all  incoming  data  with  10  millisecond  accuracy,  synchronized  to  the  timing 
pulse  of  a  GPS  receiver.  The  TEMTADS  should  be  configured  to  output  a  precise 
pulse  at  the  start  of  the  transmitter  on-time  and  an  event  identifier.  The  DAS  could 
ingest  and  time-stamp  the  pulse  and  the  event  identifier,  which  could  then  be  used  to 
assign  accurate  times  to  the  TEMTADS  data  which  would  still  be  logged  within  the 
existing  control  software  for  that  sensor. 

3)  Measurement  time  and  data-density:  The  longer  the  time-range  measured,  the  further 
the  array  will  move  during  a  single  measurement,  the  lower  the  density  of  the 
measurements  and  the  lower  the  SNR  at  the  later  time-channels.  The  TEMTADS  has 
the  ability  to  measure  the  time-decay  to  25  ms  after  transmitter  shut-off,  with  a  single 
measurement  with  that  time  range  requiring  about  100  ms  to  collect.  The  ESTCP 
MM-0504  demonstrations  at  FLBGR  and  Camp  Sibert  demonstrated  that  the  wider 
measurement  range  of  the  EM63  (to  25  ms)  resulted  in  superior  discrimination 
performance  to  the  EM61  (to  1  ms).  It’s  not  clear  as  to  whether  measurements  are 
required  all  the  way  out  to  25  ms  and  some  compromise  of  measuring  out  to,  say,  8 
ms,  may  prove  to  be  more  effective  (measurements  require  32  ms  and  hence  the  data- 
density  will  be  about  3  times  greater  than  at  25  ms  decay).  The  TEMTADS 
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measurement  characteristics  are  user-configurable  and  should  be  investigated  as  part 
of  the  proposed  follow-on  ESTCP  project.  During  that  project  we  would  also 
implement  and  investigate  the  efficacy  of  explicitly  including  the  sensor  movement 
during  measurement  in  the  polarization  tensor  inversion  process. 

4)  Signal  to  noise  ratio:  In  cued-mode,  the  data  are  stacked  to  improve  signal-to-noise 
ratio  (typically  27  repeats  of  a  negative  then  positive  transmit  cycle  are  used).  In 
dynamic  mode,  the  choice  of  the  number  of  stacks  involves  a  compromise  between 
improving  SNR  and  blurring  of  spatial  detail.  The  lower  the  number  of  stacks,  the 
less  movement  that  occurs  during  each  measurement.  The  greater  the  current  in  the 
transmitter  the  better  the  SNR,  but  the  longer  the  delay  before  the  receiver  response  to 
the  transmitter  has  decayed  sufficiently.  As  part  of  the  follow-on  project,  we 
recommend  determining  the  maximum  practical  increase  in  transmitter  current  that 
can  be  achieved  with  the  existing  TEMTADS  transmitters,  and  if  necessary,  the 
construction  of  new  transmitters  and  circuitry  to  improve  SNR  (see  6)  below). 
Movement  during  the  receive  period  also  introduces  dynamic  motion  noise  which 
reduces  the  SNR.  In  an  effort  to  combat  that  noise  source  we  recommend  including 
an  extra  reference  sensor  0.5  m  higher  than  the  existing  plane  of  sensors.  This 
reference  sensor  will  sense  approximately  the  same  ambient  noise  as  the  other 
receivers,  and  due  to  its  increased  height,  will  have  a  significantly  reduced  response 
from  any  items  of  interest. 

5)  Change  in  background  response  from  the  ground:  During  mobile  data  collection,  the 
geometry  of  the  array  relative  to  the  ground  is  constantly  changing.  Emerging  work 
from  project  SERDP  MM- 1573  indicates  that  the  resulting  variation  in  the 
background  response  can  be  significant  and  has  the  potential  to  degrade  the  accuracy 
of  any  polarization  tensors  extracted  from  the  data.  With  knowledge  of  the  array 
geometry  relative  to  the  ground,  variations  in  the  background  response  can  be 
modeled  and  removed  from  the  data.  We  recommend  testing  the  effectiveness  of 
using  an  accurate  digital  elevation  model  and  the  three-dimensional  position/ 
orientation  of  the  array  as  inputs  in  a  background  correction  procedure.  The  high- 
resolution  DEM  could  be  obtained  by  including  four  or  more  acoustic  altimeters  on 
the  front  of  the  TEMTADS  array  (a  set-up  we  have  used  successfully  with  the 
NRL/Sky  Research  HeliMag  system). 

6)  Excitation  of  multiple  axes:  To  be  able  to  discriminate  using  feature  vectors  from  all 
three  polarizations,  the  object  must  be  excited  from  a  wide-range  of  different  angles. 
In  cued-mode,  the  TEMTADS  achieves  that  by  cycling  through  each  of  the  25 
transmitters.  In  a  mobile  discrimination  mode,  there  is  not  enough  time  to  fire  each 
transmitter  separately.  In  the  simulations  section  we  investigated  the  performance  of  a 
the  TEMTADS  when  four  1  m  by  1  m  transmitters  are  incorporated  into  the  system. 
In  principal  it’s  possible  to  fire  multiple  TEMTADS  transmitters  at  the  same  time.  If 
all  five-transmitters  in  the  central  row  are  fired,  the  transmitter  would  produce  a  field 
equivalent  to  a  2  m  by  0.4  m  transmitter.  As  shown  in  Figure  34a,  this  field  would  be 
vertical  along  the  row,  and  would  be  approximately  horizontal  0.4  meters  in  front,  or 
behind  that  row  of  transmitters.  If  the  transmitter  were  fired  every,  say,  40  ms,  as  the 
array  moved  along  a  line,  subsurface  objects  would  be  excited  by  varying  amounts  in 
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the  vertical  and  the  along-line  horizontal  directions.  However,  there  would  be  very 
little  excitation  in  the  transverse  horizontal  direction.  This  is  why  systems  like  the 
MTADS  EM  array  must  conduct  two  perpendicular  passes  over  an  area  in  order  to 
excite  an  object  from  three  orthogonal  directions.  The  multiple  transmitters/receivers 
of  the  TEMTADS  along  the  direction  of  travel  open  up  the  possibility  of  exciting  the 
transverse  horizontal  components  without  the  need  to  drive  lines  in  a  perpendicular 
direction.  For  instance,  if  5  transmitters  along  a  column  (i.e.  in  the  direction  of  travel) 
are  fired,  the  field  immediately  under  that  column  would  be  vertical,  but  the  field  on 
either  side  would  be  horizontal  in  the  transverse  direction  (Figure  34a  with  rows  and 
columns  changed).  Because  there  are  5  receivers  and  transmitters  along  the  direction 
of  travel,  each  column  would  not  need  to  be  fired  as  often  as  the  central  row.  For 
instance,  a  possible  firing  sequence  would  be  to  fire  the  central  row  every-second 
measurement,  and  then  one  of  the  five  columns  for  every  other  measurement. 
Therefore,  10  measurements  would  be  required  for  one  whole  sequence  (r3,  cl,  r3, 
c2,  r3,  c3,  r3,  c4,  r3,  c5  where  r3  represents  row  number  3  etc),  and  with  a  40  ms 
measurement  time,  this  would  translate  to  0.4  seconds  or  40  cm  of  movement  at  1 
m/s.  At  that  point,  the  measurement  sequence  would  be  repeated  and  as  the  array 
continued  to  traverse  down  the  line,  objects  would  be  excited  from  the  vertical  and 
the  along-  and  across-track  directions.  Horizontal  fields  in  the  transverse  direction 
could  also  be  obtained  by  firing  one  or  more  transmitters  in  opposite  directions 
(Figure  34b).  Assuming  balanced  currents,  the  field  in  the  middle  of  the  two 
transmitters  would  be  horizontal.  While  in  principal  it’s  possible  to  fire  multiple 
TEMTADS  transmitters  at  the  same  time,  there  may  be  design  constraints  that 
preclude  that  possibility.  In  addition,  we  do  not  yet  know  if  the  transmitter  current  in 
the  smaller  40  cm  loops  will  be  large  enough  to  provide  sufficient  SNR  for  detection 
and  discrimination  at  depths  approaching  1 1  times  the  diameter  of  an  ordnance  item. 
We  recommend  conducting  some  initial  feasibility  tests  at  the  project  outset  to 
determine  if  the  existing  transmitters  can  be  utilized  in  a  one-pass  detection  and 
discrimination  mode.  If  not,  then  larger  more  powerful  transmitters  will  need  to  be 
constructed  and  aligned  in  such  a  way  that  there  is  good  excitation  in  all  three 
orthogonal  directions. 
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(a)  (b) 

Figure  34:  Cross  section  of  the  primary  field  pattern  from  the  TEMTADS  array:  (a)  when  the  center 
column  (or  row)  is  fired;  and  (b)  when  the  two  outer  columns  (or  rows)  are  fired  in  opposition  to  the 
center  column  (or  row).  The  black  dots  show  the  positions  of  the  transmitter  loops.  The  top-row 
shows  the  horizontal  and  vertical  fields  at  two  depths,  while  the  bottom  row  shows  the  field 
direction  and  strength  in  a  slice  below  the  array. 
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APPENDIX  A:  FIGURES  AND  TABLES  FROM  MONTE 
CARLO  SIMULATIONS 


In  this  appendix,  we  include  all  figures  and  tables  for  each  simulation  of  the  5  systems  under 
various  measurement  errors. 


The  top  panel  of  each  figure  displays  the  inverted  parameters  of  (k,/3, 7).  In  these  scatter 
plots,  the  kp  and  ks  represent  the  parameters  for  the  primary  and  secondary  polarizations 
that  are  chosen  as  the  maximum  and  minimum  values  among  the  three  inverted  parameters, 
and  (3P ,  3S,  yp,  7.5  are  found  accordingly.  After  the  primary  and  secondary  quantities,  the 
tertiary  ones  are  defined.  The  middle  panel  shows  the  corresponding  polarizations,  denoted 
as  the  primary  Lp(t),  secondary  Ls(t),  and  tertiary  Lt(t)  polarizations  with  true  polarization 
curves  in  shown  in  green.  The  bottom  panel  contains  three  scatter  plots,  derived  from  the 
recovered  polarizations  at  the  early  and  late  times,  for  example,  t 1  =  0.18  ms  and  f2i  =  10 
ms,  which  show  how  some  additional  information  can  be  extracted  for  discriminator  purpose. 
For  instance,  the  bottom  right  plot  measures  the  difference  between  the  secondary  and  tertiary 
polarizations  at  t\  and  f2 1-  The  small  spread  in  late  time  in  this  plot  indicates  cylindrical 
symmetry. 

The  root-mean-square  error  (RMSE)  and  the  standard  deviation  (SD)  in  the  tables  are  defined 
as 


and 


RMS  = 


\ 


TV; 


Ni  ^ 


T 7/  inlr\w.  \ 


'  i  ■  ,  \  TTltrue 

3=1 


J 


(A-2) 


to  measure  the  difference  or  the  spread  of  a  recovered  parameter  of  rrij  in  jth  simulation  from 
the  true  one  of  mtrUe-  IV*  is  the  simulation  number  for  z'th  item. 


Figures  A-l  and  A-2  and  shows  the  results  from  cued  array  when  a  =  0.06  and  0.6  in  sensor 
noise. 


Figures  A-3  and  A-4  show  the  results  from  the  towed  system  of  Single-Tx/25-Rx  in  which 
random  measurement  errors  include  sensor  noise,  perturbations  in  sensor  orientation 
with  a  standard  deviation  of  <rorien  =  2°,  as  well  as  perturbations  in  sensor  position 
with  a  standard  deviation  of  apos  =  2  cm,  5  cm. 


Figures  A-5  and  A-6  show  the  results  from  the  towed  system  of  4-Tx/25-Rx  in  which  ran¬ 
dom  measurement  errors  include  sensor  noise,  perturbations  in  sensor  orientation  with 
a  standard  deviation  of  crorien  =  2°,  as  well  as  perturbations  in  sensor  position  with  a 
standard  deviation  of  <7pos  =  2  cm  ,  5  cm. 
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Figures  A-7  and  A-8  show  the  results  from  the  cart  system  of  EM63-M  in  which  random 
measurement  errors  include  sensor  noise,  perturbations  in  sensor  orientation  with  a 
standard  deviation  of  crorieri  =  2°,  as  well  as  perturbations  in  sensor  position  with  a 
standard  deviation  of  apos  =  2  cm  ,  5  cm. 

Figures  A-9  and  A- 10  show  the  results  from  the  cart  system  of  EM63  in  which  random 
measurement  errors  include  sensor  noise,  perturbations  in  sensor  orientation  with  a 
standard  deviation  of  o-0nen  =  2°,  as  well  as  perturbations  in  sensor  position  with  a 
standard  deviation  of  <7pos  =  2  cm  ,  5  cm. 

Figures  A- 11  and  A- 13  show  the  results  from  the  cart  system  of  EM63  in  which  random 
measurement  errors  are  only  from  perturbations  in  sensor  position  with  a  standard  de¬ 
viation  of  u pos  =  1  cm,  2  cm,  and  5  cm. 

Table  A-l  list  the  statistical  measures  of  Mean,  SD,  Min(imum),  Max(imum),  and  RMSE 
for  each  set  of  recovered  parameters  of  ( k ,  f3, 7)  in  the  three  polarizations  using  the 
cued  array  with  different  levels  of  sensor  noise. 

Table  A- 2  list  the  statistical  measures  of  Mean,  SD,  Min(imum),  Max(imum),  and  RMSE  for 
each  set  of  recovered  parameters  of  (k,  /?,  7)  in  the  three  polarizations  using  the  Single- 
Tx/25-Rx  system  with  different  levels  of  measurement  errors  arising  from  sensor  noise 
as  well  as  perturbations  in  sensor  position  and  orientation. 

Table  A-3  list  the  statistical  measures  of  Mean,  SD,  Min(imum),  Max(imum),  and  RMSE 
for  each  set  of  recovered  parameters  of  (k,  (3, 7)  in  the  three  polarizations  using  Four- 
Tx/25-Rx  system  with  different  levels  of  measurement  errors  arising  from  sensor  noise 
as  well  as  perturbations  in  sensor  position  and  orientation. 

Table  A-4  list  the  statistical  measures  of  Mean,  SD,  Min(imum),  Max(imum),  and  RMSE 
for  each  set  of  recovered  parameters  of  (k,  f3, 7)  in  the  three  polarizations  using  the 
EM63-M  system  with  different  levels  of  measurement  errors  arising  from  sensor  noise 
as  well  as  perturbations  in  sensor  position  and  orientation. 

Table  A-5  list  the  statistical  measures  of  Mean,  SD,  Min(imum,  Max(imum),  and  RMSE  for 
each  set  of  recovered  parameters  of  (k,  (3, 7)  in  the  three  polarizations  using  the  EM63 
system  with  different  levels  of  measurement  errors  arising  from  sensor  noise  as  well  as 
perturbations  in  sensor  position  and  orientation. 

Table  A-6  list  the  statistical  measures  of  Mean,  SD,  Min(imum),  Max(imum),  and  RMSE 
for  each  set  of  recovered  parameters  of  (k,  (3 , 7)  in  the  three  polarizations  using  Single- 
Tx/25-Rx  system  with  measurement  errors  solely  arising  from  different  levels  of  per¬ 
turbations  in  sensor  position. 
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Figure  A-l:  The  results  of  three-polarization 
sensor  noise  of  ^  *  randn. 
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Figure  A-2:  The  results  of  three-polarization 
sensor  noise  of  *  randn. 
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Figure  A-3:  The  results  of  three-polarization  inversion  of  noisy  data  using  Single-Tx/25-Rx 
system.  Measurement  errors  arising  from  sensor  noise  *  randn)  and  perturbations  in 
sensor  position  and  orientation  with  crorien  =  2°,  <jpos  =  2  cm. 
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Figure  A-4:  The  results  of  three-polarization  inversion  of  noisy  data  using  Single-Tx/25-Rx 
system.  Measurement  errors  arising  from  sensor  noise  (^|  *  randn)  and  perturbations  in 

sensor  position  and  orientation  with  crorien  =  2°,  apos  =  5  cm. 
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Figure  A-5:  The  results  of  three-polarization  inversion  of  noisy  data  using  4-Tx/25-Rx  sys¬ 
tem.  Measurement  errors  arising  from  sensor  noise  (^|  *  randn)  and  perturbations  in  sensor 

position  and  orientation  with  rronen  =2°,  <7pOS  =  2  cm. 
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Figure  A-6:  The  results  of  three-polarization  inversion  of  noisy  data  using  4-Tx/25-Rx  sys¬ 
tem.  Measurement  errors  arising  from  sensor  noise  (^|  *  randn)  and  perturbations  in  sensor 

position  and  orientation  with  rroneri  =2°,  <7pOS  =  5  cm. 
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Figure  A-7:  The  results  of  three-polarization  inversion  of  noisy  data  using  EM63-M  sys¬ 
tem.  Measurement  errors  arising  from  sensor  noise  *  randn)  and  perturbations  in  sensor 
position  and  orientation  with  rronen  =2°,  <7pos  =  2  cm. 
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Figure  A-8:  The  results  of  three-polarization  inversion  of  noisy  data  using  EM63-M  sys¬ 
tem.  Measurement  errors  arising  from  sensor  noise  (^|  *  randn)  and  perturbations  in  sensor 

position  and  orientation  with  rronen  =2°,  <7pOS  =  5  cm. 
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Figure  A-9:  The  results  of  three-polarization  inversion  of  noisy  data  using  EM63  system. 
Measurement  errors  arising  from  sensor  noise  (77=  *  randn)  and  perturbations  in  sensor  posi¬ 
tion  and  orientation  with  crorien  =  2°,  crpos  =  2  cm. 
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Figure  A- 10:  The  results  of  three-polarization  inversion  of  noisy  data  using  EM63  system. 
Measurement  errors  arising  from  sensor  noise  (77=  *  randn)  and  perturbations  in  sensor  posi¬ 
tion  and  orientation  with  crorien  =  2°,  crpos  =  5  cm. 
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Figure  A-ll:  The  results  of  three-polarization  inversion  of  noisy  data  using  EM63  system. 
Measurement  errors  only  arising  from  perturbations  in  sensor  position  with  <jpos  =  1  cm. 


A-13 


September  2008 


MM-1571:  Next  Generation  Data  Collection  Platform 


Sky  Research  Inc. 


O  37  mm 
0  81  mm 

□  4.2  in 


O  37  mm 
<>  81  mm 

□  4.2  in 


MV 


MV 


O  37  mm 
0  81  mm 

2  □  4.2  in 


^  0- 


10  10  10 

lL,<'2,>-MWI/Ls,,rue<W 


10 


Figure  A- 12:  The  results  of  three-polarization  inversion  of  noisy  data  using  EM63  system. 
Measurement  errors  only  arising  from  perturbations  in  sensor  position  with  <jpos  =  2  cm. 
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Figure  A- 13:  The  results  of  three-polarization  inversion  of  noisy  data  using  EM63  system. 
Measurement  errors  only  arising  from  perturbations  in  sensor  position  with  <jpos  =  5  cm. 
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Table  A-l:  Statistical  results  for  3 -polarization  inversion  of  noisy  data  using  cued  array. 


Parameters 

kp 

ks 

h 

Pp 

Ps 

Pt 

7 P 

7s 

7 1 

Sensor  noise:  0.06  jyft. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.83 

0.62 

0.62 

0.44 

1.16 

1.15 

4.25 

5.79 

5.66 

SD 

0.06 

0.01 

0.01 

0.00 

0.00 

0.01 

0.01 

0.33 

0.25 

Min 

4.52 

0.59 

0.60 

0.44 

1.15 

1.13 

4.21 

5.15 

5.26 

Max 

4.85 

0.63 

0.63 

0.45 

1.17 

1.17 

4.28 

7.25 

6.85 

RMSE 

0.01 

0.01 

0.01 

0.00 

0.00 

0.00 

0.00 

0.06 

0.04 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

23.88 

9.71 

9.77 

0.62 

0.79 

0.78 

3.73 

1.88 

1.86 

SD 

0.04 

0.05 

0.07 

0.00 

0.00 

0.00 

0.02 

0.02 

0.02 

Min 

23.80 

9.55 

9.72 

0.62 

0.78 

0.77 

3.65 

1.85 

1.79 

Max 

24.04 

9.78 

9.91 

0.62 

0.80 

0.79 

3.78 

1.94 

1.89 

RMSE 

0.00 

0.01 

0.01 

0.00 

0.00 

0.00 

0.01 

0.01 

0.01 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

98.25 

22.53 

22.66 

0.53 

0.95 

0.95 

4.33 

3.84 

3.77 

SD 

0.22 

0.25 

0.14 

0.00 

0.01 

0.01 

0.02 

0.16 

0.07 

Min 

97.18 

21.22 

22.54 

0.53 

0.95 

0.92 

4.23 

3.77 

3.45 

Max 

98.72 

22.63 

23.23 

0.54 

0.98 

0.95 

4.40 

4.69 

3.83 

RMSE 

0.00 

0.01 

0.01 

0.00 

0.01 

0.01 

0.01 

0.04 

0.02 

Sensor  noise:  0.6  j\[t. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.76 

0.60 

0.64 

0.44 

1.16 

1.14 

4.29 

6.48 

5.98 

SD 

0.26 

0.03 

0.06 

0.02 

0.02 

0.04 

0.28 

2.04 

3.26 

Min 

3.37 

0.54 

0.60 

0.42 

1.08 

0.95 

3.98 

3.60 

3.40 

Max 

4.85 

0.63 

0.94 

0.52 

1.19 

1.17 

5.83 

12.92 

23.78 

RMSE 

0.05 

0.04 

0.09 

0.03 

0.02 

0.03 

0.07 

0.36 

0.57 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

23.51 

9.37 

9.86 

0.62 

0.79 

0.77 

3.66 

1.96 

1.79 

SD 

0.80 

0.59 

0.48 

0.01 

0.02 

0.03 

0.26 

0.22 

0.22 

Min 

20.42 

8.19 

8.62 

0.58 

0.75 

0.66 

2.91 

1.66 

1.31 

Max 

24.31 

9.98 

11.45 

0.64 

0.85 

0.81 

4.33 

2.51 

2.57 

RMSE 

0.03 

0.06 

0.05 

0.02 

0.03 

0.04 

0.07 

0.12 

0.12 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

97.33 

22.17 

22.81 

0.53 

0.96 

0.94 

4.33 

3.91 

3.61 

SD 

2.05 

0.68 

0.56 

0.01 

0.01 

0.02 

0.17 

0.27 

0.32 

Min 

92.43 

20.88 

21.28 

0.50 

0.94 

0.90 

3.66 

3.58 

2.79 

Max 

100.98 

22.63 

24.08 

0.54 

0.99 

0.96 

4.95 

4.70 

3.92 

RMSE 

0.02 

0.03 

0.02 

0.01 

0.01 

0.02 

0.04 

0.07 

0.08 
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Table  A-2:  Statistical  results  for  3 -polarization  inversion  of  noisy  data  using  Single-Tx/25- 
Rx. 


Parameters 

kp 

ks 

kt 

Pp 

Ps 

A 

7 P 

7 « 

7 1 

Sensor  noise  of  ^4 

*  randn  and  sensor  errors  with  crorien 

=  2U, 

crp0s  =  2  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.69 

0.60 

0.72 

0.43 

1.15 

1.07 

4.59 

7.83 

8.17 

SD 

0.28 

0.09 

0.18 

0.03 

0.07 

0.13 

2.37 

6.20 

8.92 

Min 

3.97 

0.44 

0.57 

0.33 

0.97 

0.77 

2.18 

1.37 

1.00 

Max 

5.45 

0.84 

1.32 

0.49 

1.31 

1.20 

17.42 

24.21 

39.04 

RMSE 

0.06 

0.14 

0.29 

0.06 

0.06 

0.11 

0.56 

1.09 

1.57 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

22.72 

8.72 

9.97 

0.62 

0.81 

0.75 

4.09 

3.62 

1.89 

SD 

1.95 

1.42 

1.27 

0.03 

0.06 

0.07 

1.60 

4.84 

4.84 

Min 

18.33 

5.98 

7.28 

0.50 

0.69 

0.58 

2.03 

1.13 

1.00 

Max 

28.01 

10.43 

12.87 

0.70 

0.97 

0.87 

11.59 

18.16 

4.84 

RMSE 

0.08 

0.15 

0.13 

0.05 

0.08 

0.09 

0.43 

2.59 

2.59 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

93.46 

21.01 

22.89 

0.53 

0.96 

0.94 

4.38 

4.61 

4.00 

SD 

8.53 

2.09 

1.61 

0.03 

0.03 

0.03 

0.55 

2.69 

0.83 

Min 

64.86 

16.06 

18.73 

0.44 

0.90 

0.83 

2.48 

3.21 

2.03 

Max 

105.63 

23.11 

27.72 

0.60 

1.04 

1.05 

6.70 

17.41 

6.34 

RMSE 

0.09 

0.09 

0.07 

0.05 

0.03 

0.04 

0.13 

0.71 

0.22 

Sensor  noise  of  ^4 

*  randn  and  sensor  errors  with  crorien 

=  2U, 

^pos  —  5  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.50 

0.57 

0.82 

0.42 

1.15 

1.01 

5.09 

6.18 

10.81 

SD 

0.51 

0.19 

0.31 

0.03 

0.10 

0.17 

2.52 

5.54 

10.20 

Min 

3.99 

0.32 

0.61 

0.36 

0.85 

0.69 

3.98 

1.00 

1.18 

Max 

6.08 

1.44 

1.69 

0.47 

1.40 

1.13 

17.28 

24.87 

31.51 

RMSE 

0.11 

0.30 

0.50 

0.07 

0.09 

0.15 

0.59 

0.98 

1.80 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

22.46 

8.15 

9.82 

0.61 

0.81 

0.77 

4.11 

3.64 

3.70 

SD 

2.71 

1.88 

1.69 

0.03 

0.07 

0.09 

2.83 

4.86 

4.86 

Min 

17.17 

5.55 

6.44 

0.52 

0.72 

0.61 

2.24 

1.38 

1.02 

Max 

28.14 

10.30 

12.93 

0.71 

0.99 

0.97 

19.75 

22.72 

24.91 

RMSE 

0.11 

0.19 

0.17 

0.05 

0.08 

0.11 

0.76 

2.60 

2.60 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

90.37 

19.65 

22.73 

0.51 

0.95 

0.92 

4.35 

4.89 

4.43 

SD 

14.14 

3.72 

2.18 

0.03 

0.05 

0.06 

0.43 

2.42 

1.45 

Min 

56.56 

12.95 

15.16 

0.43 

0.82 

0.75 

3.42 

2.66 

2.24 

Max 

104.26 

22.98 

26.99 

0.56 

1.01 

1.02 

6.26 

11.26 

8.40 

RMSE 

0.14 

0.16 

0.10 

@U©6 

0.05 

0.06 

0.10 

0.64 
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Table  A-3:  Statistical  results  for  3 -polarization  inversion  of  noisy  data  using  4-Tx/25-Rx. 


Parameters 

kp 

ks 

kt 

Pp 

Ps 

Pt 

7 P 

7s 

7 1 

Sensor  noise  of 

*  randn  and  sensor  errors  with  crorien 

=  2U, 

crpos  =  2  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.64 

0.58 

0.66 

0.43 

1.16 

1.11 

4.36 

7.72 

5.63 

SD 

0.26 

0.06 

0.08 

0.01 

0.05 

0.06 

0.38 

5.08 

2.73 

Min 

4.29 

0.46 

0.56 

0.41 

1.10 

0.97 

3.93 

2.32 

1.57 

Max 

4.98 

0.69 

0.94 

0.49 

1.30 

1.17 

6.29 

21.78 

13.10 

RMSE 

0.05 

0.10 

0.12 

0.03 

0.04 

0.05 

0.09 

0.90 

0.48 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

22.81 

9.06 

10.00 

0.62 

0.79 

0.76 

3.91 

2.44 

1.79 

SD 

1.54 

1.17 

0.83 

0.01 

0.06 

0.05 

0.65 

2.20 

2.20 

Min 

19.55 

5.72 

7.44 

0.58 

0.68 

0.66 

3.44 

1.23 

1.26 

Max 

24.21 

10.70 

11.75 

0.66 

0.98 

0.87 

7.15 

11.68 

3.15 

RMSE 

0.06 

0.12 

0.08 

0.02 

0.07 

0.06 

0.17 

1.18 

1.18 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

95.08 

20.82 

22.80 

0.52 

0.97 

0.93 

4.33 

5.13 

3.72 

SD 

5.41 

2.45 

1.20 

0.02 

0.04 

0.04 

0.43 

3.82 

0.56 

Min 

83.69 

16.48 

19.78 

0.47 

0.93 

0.81 

2.88 

3.46 

2.26 

Max 

101.89 

23.11 

24.60 

0.56 

1.09 

0.98 

6.01 

24.03 

5.09 

RMSE 

0.05 

0.11 

0.05 

0.04 

0.04 

0.04 

0.10 

1.01 

0.15 

Sensor  noise  of  ^4 

*  randn  and  sensor  errors  with  <7orien 

=  2U, 

^pos  5  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.28 

0.56 

0.71 

0.43 

1.12 

1.05 

4.65 

8.97 

9.43 

SD 

0.73 

0.12 

0.14 

0.02 

0.08 

0.12 

0.64 

8.67 

9.87 

Min 

3.22 

0.31 

0.48 

0.38 

0.95 

0.90 

3.78 

1.25 

1.17 

Max 

5.03 

0.81 

1.07 

0.46 

1.29 

1.24 

6.37 

30.13 

43.16 

RMSE 

0.15 

0.19 

0.23 

0.05 

0.07 

0.11 

0.15 

1.53 

1.74 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

22.40 

8.70 

10.13 

0.61 

0.79 

0.75 

4.12 

2.65 

2.10 

SD 

2.52 

1.62 

1.54 

0.03 

0.08 

0.08 

1.50 

1.74 

1.74 

Min 

18.15 

6.46 

7.02 

0.48 

0.62 

0.52 

3.24 

1.54 

1.10 

Max 

26.19 

12.22 

14.27 

0.70 

0.96 

0.91 

12.08 

8.77 

4.89 

RMSE 

0.11 

0.17 

0.16 

0.06 

0.10 

0.10 

0.40 

0.93 

0.93 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

88.54 

19.85 

22.80 

0.51 

0.95 

0.90 

4.52 

5.00 

3.73 

SD 

13.77 

4.06 

2.16 

0.05 

0.10 

0.07 

0.35 

3.35 

0.86 

Min 

61.28 

8.65 

16.45 

0.30 

0.56 

0.71 

3.61 

2.52 

1.22 

Max 

100.85 

23.19 

27.76 

0.58 

1.24 

0.98 

5.20 

19.38 

5.52 

RMSE 

0.14 

0.18 

0.10 

0.09 

0.11 

0.08 

0.08 

0.88 

0.23 
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Table  A-4:  Statistical  results  for  3 -polarization  inversion  of  noisy  data  using  EM63-M. 


Parameters 

kp 

ks 

kt 

Pp 

Pa 

Pt 

7 P 

Is 

It 

Sensor  noise  of 

*  randn  and  sensor  errors  with  crorien 

=  2U, 

crpos  =  2  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.97 

0.59 

0.79 

0.43 

1.12 

1.06 

4.39 

5.58 

6.44 

SD 

0.58 

0.12 

0.28 

0.04 

0.15 

0.16 

0.53 

4.62 

6.23 

Min 

4.12 

0.29 

0.42 

0.32 

0.68 

0.66 

2.71 

1.00 

1.00 

Max 

6.09 

0.79 

1.69 

0.52 

1.50 

1.37 

5.98 

23.84 

30.01 

RMSE 

0.12 

0.19 

0.44 

0.09 

0.13 

0.14 

0.12 

0.81 

1.10 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

24.48 

8.50 

10.36 

0.62 

0.82 

0.78 

5.74 

2.70 

2.62 

SD 

4.45 

1.80 

2.04 

0.05 

0.09 

0.08 

8.38 

1.71 

1.71 

Min 

15.09 

6.10 

6.85 

0.49 

0.71 

0.64 

1.68 

1.32 

1.22 

Max 

41.25 

11.19 

16.63 

0.74 

1.01 

1.01 

48.76 

8.90 

16.84 

RMSE 

0.19 

0.19 

0.21 

0.08 

0.11 

0.10 

2.24 

0.92 

0.92 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

97.48 

21.21 

25.19 

0.53 

0.97 

0.90 

4.99 

5.47 

3.39 

SD 

8.30 

3.29 

3.48 

0.05 

0.07 

0.08 

1.76 

4.09 

0.88 

Min 

81.64 

11.91 

20.18 

0.30 

0.73 

0.75 

3.49 

2.01 

1.81 

Max 

110.68 

26.75 

31.28 

0.64 

1.08 

1.05 

11.74 

17.45 

4.46 

RMSE 

0.08 

0.15 

0.15 

0.10 

0.07 

0.09 

0.41 

1.08 

0.23 

Sensor  noise  of  ^4 

*  randn  and  sensor  errors  with  <7orien 

=  2U, 

^pos  5  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

5.01 

0.52 

0.89 

0.43 

1.08 

0.96 

5.21 

8.03 

9.60 

SD 

1.79 

0.22 

0.41 

0.07 

0.22 

0.28 

1.89 

9.87 

11.35 

Min 

1.56 

0.21 

0.52 

0.30 

0.45 

0.33 

3.39 

1.00 

1.00 

Max 

12.29 

1.29 

1.81 

0.64 

1.50 

1.30 

11.69 

43.68 

47.08 

RMSE 

0.37 

0.36 

0.65 

0.17 

0.19 

0.24 

0.44 

1.74 

2.00 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

26.51 

8.59 

12.30 

0.61 

0.81 

0.74 

4.84 

4.07 

2.61 

SD 

7.28 

2.42 

6.06 

0.04 

0.14 

0.12 

3.30 

6.22 

6.22 

Min 

12.75 

4.79 

6.06 

0.53 

0.53 

0.58 

2.56 

1.00 

1.00 

Max 

43.18 

14.23 

32.51 

0.71 

1.25 

0.96 

18.09 

27.15 

12.07 

RMSE 

0.31 

0.25 

0.62 

0.06 

0.18 

0.15 

0.88 

3.32 

3.32 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

96.71 

20.05 

25.13 

0.53 

0.96 

0.90 

5.49 

6.07 

3.53 

SD 

23.91 

6.20 

6.03 

0.06 

0.10 

0.12 

3.43 

8.49 

1.78 

Min 

32.51 

5.03 

14.44 

0.35 

0.79 

0.59 

3.54 

1.50 

1.00 

Max 

144.89 

30.87 

37.86 

0.72 

1.31 

1.13 

20.83 

50.00 

10.28 

RMSE 

0.24 

0.27 

0.27 

0.12 

0.10 

0.12 

0.79 

2.23 

0.47 
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Table  A-5:  Statistical  results  for  3 -polarization  inversion  of  noisy  data  using  EM63. 


Parameters 

kp 

ks 

kt 

Pp 

Pa 

Pt 

7 P 

Is 

It 

Sensor  noise  of 

*  randn  and  sensor  errors  with  crorien 

=  2U, 

crpos  =  2  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.69 

0.57 

0.79 

0.44 

1.13 

1.06 

4.94 

6.25 

9.15 

SD 

1.45 

0.16 

0.29 

0.05 

0.14 

0.18 

3.49 

4.73 

8.76 

Min 

0.66 

0.32 

0.41 

0.30 

0.75 

0.61 

3.24 

1.00 

1.00 

Max 

8.38 

1.03 

1.43 

0.62 

1.41 

1.44 

24.11 

22.09 

26.85 

RMSE 

0.30 

0.26 

0.46 

0.11 

0.12 

0.16 

0.82 

0.83 

1.54 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

29.60 

7.18 

13.07 

0.62 

0.84 

0.69 

5.03 

4.52 

3.82 

SD 

21.48 

3.59 

13.15 

0.07 

0.20 

0.15 

8.11 

8.71 

8.71 

Min 

9.86 

0.43 

5.86 

0.42 

0.30 

0.30 

1.12 

1.00 

1.00 

Max 

107.31 

11.62 

82.52 

0.75 

1.50 

0.83 

50.00 

48.75 

49.99 

RMSE 

0.90 

0.37 

1.35 

0.11 

0.25 

0.20 

2.17 

4.65 

4.65 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

97.63 

20.89 

27.05 

0.53 

0.97 

0.89 

4.55 

6.49 

4.38 

SD 

31.76 

6.86 

12.32 

0.04 

0.09 

0.12 

1.35 

7.35 

5.38 

Min 

58.71 

8.65 

10.45 

0.44 

0.80 

0.59 

2.22 

1.40 

1.00 

Max 

228.19 

50.93 

71.36 

0.63 

1.23 

1.04 

10.56 

37.34 

33.74 

RMSE 

0.32 

0.30 

0.55 

0.07 

0.10 

0.13 

0.31 

1.93 

1.42 

Sensor  noise  of  ^4 

*  randn  and  sensor  errors  with  <7orien 

=  2U, 

^pos  5  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.02 

0.47 

0.83 

0.45 

1.04 

0.92 

5.30 

8.65 

8.82 

SD 

2.64 

0.31 

0.36 

0.08 

0.25 

0.29 

2.69 

9.58 

7.62 

Min 

0.61 

0.14 

0.30 

0.31 

0.47 

0.41 

3.03 

1.69 

1.00 

Max 

10.92 

1.53 

1.87 

0.68 

1.50 

1.17 

17.02 

44.72 

32.22 

RMSE 

0.55 

0.50 

0.58 

0.17 

0.22 

0.25 

0.63 

1.69 

1.34 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

26.45 

5.91 

10.82 

0.62 

0.76 

0.72 

4.86 

3.31 

3.65 

SD 

34.96 

4.69 

9.74 

0.07 

0.13 

0.13 

3.78 

4.39 

4.39 

Min 

4.00 

1.18 

3.31 

0.45 

0.52 

0.30 

1.19 

1.00 

1.00 

Max 

205.94 

11.50 

56.73 

0.88 

1.11 

0.87 

22.55 

25.50 

29.22 

RMSE 

1.46 

0.48 

1.00 

0.11 

0.17 

0.16 

1.01 

2.35 

2.35 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

95.52 

16.03 

25.51 

0.55 

0.98 

0.88 

7.68 

8.38 

4.24 

SD 

44.10 

8.03 

10.35 

0.06 

0.15 

0.12 

10.44 

11.44 

1.56 

Min 

12.85 

8.00 

8.68 

0.48 

0.47 

0.60 

2.96 

1.06 

1.16 

Max 

234.68 

25.17 

49.55 

0.72 

1.34 

1.06 

49.73 

44.09 

7.71 

RMSE 

0.45 

0.36 

0.46 

0.11 

0.15 

0.12 

2.41 

3.01 

0.41 
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Table  A-6:  Statistical  results  for  3 -polarization  inversion  using  EM63.  Only  sensor  positional 
errors. 


Parameters 

kp 

ks 

kt 

Pp 

Pa 

A 

7 P 

7s 

7 1 

Perturbations  in  sensor  position  with  crpos  = 

=  1  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.91 

0.58 

0.69 

0.44 

1.17 

1.13 

4.28 

5.07 

6.40 

SD 

0.59 

0.07 

0.08 

0.01 

0.03 

0.04 

0.09 

1.96 

2.78 

Min 

3.54 

0.48 

0.58 

0.42 

1.11 

1.03 

4.11 

1.29 

1.50 

Max 

6.61 

0.74 

0.84 

0.46 

1.21 

1.18 

4.55 

8.83 

14.38 

RMSE 

0.12 

0.11 

0.14 

0.02 

0.02 

0.04 

0.02 

0.35 

0.49 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

24.31 

9.36 

10.21 

0.62 

0.79 

0.78 

3.69 

1.92 

1.93 

SD 

3.95 

0.99 

1.27 

0.01 

0.01 

0.01 

0.10 

0.22 

0.22 

Min 

9.44 

6.69 

6.84 

0.61 

0.76 

0.75 

3.41 

1.49 

1.17 

Max 

31.11 

10.91 

12.81 

0.63 

0.82 

0.81 

3.79 

2.64 

2.78 

RMSE 

0.17 

0.10 

0.13 

0.01 

0.02 

0.01 

0.03 

0.12 

0.12 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

96.01 

21.32 

22.97 

0.53 

0.95 

0.94 

4.33 

3.81 

3.91 

SD 

8.83 

2.16 

1.46 

0.01 

0.01 

0.01 

0.05 

0.24 

0.40 

Min 

78.59 

16.39 

19.63 

0.52 

0.92 

0.91 

4.25 

3.23 

2.82 

Max 

118.28 

24.79 

26.14 

0.56 

0.98 

0.97 

4.53 

4.28 

5.13 

RMSE 

0.09 

0.10 

0.06 

0.01 

0.01 

0.01 

0.01 

0.06 

0.11 

Perturbations  in  sensor  position  with  crpos  = 

=  2  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

4.54 

0.54 

0.73 

0.43 

1.16 

1.08 

4.32 

5.81 

10.30 

SD 

1.41 

0.16 

0.22 

0.02 

0.06 

0.10 

0.19 

4.28 

7.85 

Min 

1.53 

0.35 

0.51 

0.38 

0.94 

0.87 

3.93 

1.02 

1.00 

Max 

9.62 

1.12 

1.57 

0.47 

1.28 

1.18 

4.80 

17.31 

29.36 

RMSE 

0.29 

0.25 

0.35 

0.04 

0.05 

0.09 

0.04 

0.75 

1.38 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

21.53 

8.14 

9.67 

0.62 

0.80 

0.77 

4.24 

2.38 

2.14 

SD 

7.22 

2.44 

2.28 

0.02 

0.03 

0.06 

3.13 

1.44 

1.44 

Min 

5.25 

3.97 

4.23 

0.60 

0.77 

0.53 

3.43 

1.58 

1.40 

Max 

34.39 

11.47 

16.24 

0.71 

0.91 

0.88 

21.71 

7.92 

5.24 

RMSE 

0.30 

0.25 

0.23 

0.03 

0.04 

0.07 

0.84 

0.77 

0.77 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

105.60 

22.57 

26.15 

0.53 

0.95 

0.93 

4.38 

3.69 

3.88 

SD 

23.94 

3.94 

7.85 

0.01 

0.03 

0.03 

0.12 

0.73 

0.83 

Min 

72.64 

18.19 

20.64 

0.51 

0.89 

0.80 

4.19 

1.98 

1.44 

Max 

194.82 

40.14 

62.87 

0.56 

1.02 

0.97 

4.69 

5.38 

6.17 

RMSE 

0.24 

0.17 

0.35 

0.02 

0.03 

0.04 

0.03 

0.19 

0.22 
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Cont.  of  Table  A-6. 


Parameters 

kp 

ks 

kt 

Pp 

Ps 

Pt 

7 P 

7s 

It 

Perturbations  in  sensor  position  with  crpos  = 

=  5  cm. 

37  mm 

True 

4.84 

0.62 

0.62 

0.44 

1.15 

1.15 

4.25 

5.67 

5.67 

Mean 

3.53 

0.44 

0.78 

0.45 

1.01 

0.95 

5.87 

7.53 

12.28 

SD 

2.43 

0.33 

0.42 

0.06 

0.29 

0.26 

4.27 

6.33 

12.74 

Min 

0.43 

0.15 

0.30 

0.37 

0.34 

0.40 

3.98 

1.00 

1.04 

Max 

7.43 

1.73 

2.14 

0.70 

1.46 

1.18 

23.87 

22.51 

48.89 

RMSE 

0.50 

0.53 

0.67 

0.14 

0.25 

0.23 

1.00 

1.12 

2.25 

81  mm 

True 

23.88 

9.73 

9.73 

0.62 

0.79 

0.79 

3.73 

1.87 

1.87 

Mean 

30.09 

6.49 

10.32 

0.63 

0.78 

0.75 

3.69 

2.84 

2.98 

SD 

31.56 

4.65 

5.27 

0.03 

0.14 

0.05 

0.56 

2.09 

2.09 

Min 

3.88 

1.44 

3.24 

0.58 

0.55 

0.67 

2.35 

1.00 

1.45 

Max 

134.34 

13.57 

24.40 

0.70 

1.27 

0.86 

4.74 

10.45 

14.43 

RMSE 

1.32 

0.48 

0.54 

0.05 

0.17 

0.06 

0.15 

1.12 

1.12 

4.2  in 

True 

98.30 

22.60 

22.60 

0.53 

0.95 

0.95 

4.33 

3.80 

3.80 

Mean 

133.23 

19.50 

32.38 

0.55 

0.95 

0.86 

4.80 

4.57 

5.38 

SD 

142.04 

8.91 

23.03 

0.09 

0.14 

0.13 

0.80 

2.63 

5.27 

Min 

19.19 

6.26 

11.39 

0.44 

0.64 

0.57 

3.83 

1.00 

1.35 

Max 

679.48 

45.04 

98.00 

0.87 

1.49 

0.97 

6.60 

13.93 

24.81 

RMSE 

1.45 

0.39 

1.02 

0.17 

0.14 

0.13 

0.18 

0.69 

1.39 
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